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PREFACE 

■"PHE present worb is intended for the use of students who 
are learning the practice of Interpolation and Numerical 
Integration. The advantages of a practical kno wedge of this 
part of Mathematics are so obvious that it is needless to insist 
on them here ; and these subjects form an important part of the 
course in the modern Mathematical Laboratory. There are, 
however, so many claims on the time of students that the extent 
of this course, as of all others, must be kept within narrow 
limits : and it has therefore been necessary to restrict the treat- 
ment to the most central and indispensable theorems. A large 
number of numerical illustrations and examples has been given 
of a kind likely to occur in the applications of Mathematics. 

My thanks are due to Professor Whittaker and to my 
colleague, Mr E. M. Horsburgh, M.A., B.Sc, Assoc. M.Inst. C.E., 
for their valuable criticisms and suggestions. 

D. G. 

The Mathematical Laboratory, 

University of Edinburgh, 

Jidy 1915. 
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PRELIMINARY NOTE ON COMPUTATION. 

DEFORE introducing the formulae required in interpolation 
and numerical integration, it may be advisable to mention 
a few points ■which will facilitate the work of computing. 

Where computation is performed to any considerable extent, 
computer's desks will be found useful. Those used in the mathe- 
matical laboratory of the University of Edinburgh are 3' 0" wide, 
r 9" from front to back, and 2' 6|" high. They contain a locker, 
in which computing paper can be kept without being foMed, and 
a cupboard for books, and are fitted with a strong adjustable 
book-rest. Thus the computer can command a large space and 
utilize it for books, papers, drawing-board, arithmometer, or 
instruments. Each desk is supplied with a copy of Barlow's tables 
(which give the square, square root, cube, cube root, and the 
reciprocal of all numbers up to 10,000), a copy of Crelle's multi- 
plication table (which gives at sight the product of any two 
numbers each less than 1000), and with tables giving the values 
of the trigonometric functions and logarithms. These may, 
of course, be supplemented by a slide rule, or any of the various 
calculating machines* now in use, and such books of tables as 
bear particularly on the subject in question. 

Success in computation depends partly on the proper choice of 
a formula and partly on a neat and methodical arrangement of the 
work. For the latter computing paper is essential, A convenient 
size for such a paper is 26" by 16"; this should be divided 
by faint ruling into \" squares, each of which is capable of 



* Tor descriptions of these see Modern Instruments and Methods of Cal- 
culation, edited by E. M. Horsburgh. London : 6. Bell & Sons. 1914. 
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holding two digits. It will be found conducive towards accuracy 
and speed if, instead of taking down a number one digit at <a time, 
the computer takes it down two digits at a time, e.g., instead of 
taking down the individual digits 2, 0, 4, 7, 6, 3, it will be found 
better to group them together, as 20, 47, 63. Every computation 
should be performed with ink in preference to pencil; this not 
only ensures a much more lasting record of the work but also 
prevents eye-strain and fatigue. 

It need hardly be emphasized that seven-place accuracy cannot 
be obtained by the use of four-figure tables. A fairly safe rule is 
to make use of tables containing one digit more than the accuracy 
required ; to use more accurate tables than these is simply to 
increase the amount of labour without increasing the efficiency. 
For the same reason contracted methods of multiplication and 
division should be adopted when machines and tables are not at 
hand. 



CHAPTER I. 

SOME THEOREMS IN THE CALCULUS OF 
FINITE DIFFERENCES. 

1. Differences. 

Let ^/{'x) denote the increment of a function f(x) correspond- 
ing to a given constant increment w of the variable x, so that 

^f{x)=f{x + w)-f{x). 

The expression A/{x), which is usually called the First Differ- 
ence of the function f{x) corresponding to the constant increment 
w of the variable x, will, in general, be another function of x, and 
as such will have a first difference which will be obtained by 
operating on A_/(a;) with A. Calling the result of this operation 
A'/(a;), we have 

A^f{x) = A{Af{x)} 

= A{f{x + w}-/{x)} 

= {f{x+2w) -/{x + w)} - {/{x + w)-f{x)} 

=f{x + 2w) - 2/(03 + w) +/{x). 

This is termed the Second Difference of f{x) corresponding to the 
constant increment w of x. By repeating the process we obtain 
in turn AV(a;), AY(a!), ... A"/{x), which are called the 3rd, 
4th, ... wth differences of /{x) corresponding to this constant 
increment w. 

Note. — It is always possible by means of a suitable trans- 
formation to make the value of w unity. 

For since iv is the increment of x in f{x), then unity is the 

increment of — . Replacing x in f{x) by wy, we obtain a new 
w 

function F (y), which is such that 

AF{y) = F{y+l)-F{y). 
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Example 1.— Tabulate and difference the suooessive values of f(x) = ^ 
from a;=110 to a;=U8. 



X 


X' 


A 


A" 


110 


1331000 


36631 




HI 


1367631 


37297 


666 


112 


1404928 


37969 


672 


113 


1442897 


38647 


678 


114 


1481544 


39331 


684 


115 


1520875 


40021 


690 


116 


1560896 


40717 


696 


117 


1601613 


41419 


702 


118 


1643032 







6? 

6 
6 
6 
6 
6 
6 



A* 












In the above table the column headed x^ was obtained from Barlow's 
Tables. The column A gives the first differences, being obtained by sub- 
tracting each entry of the column a;' from the entry following it, e.g. if 
a; = 114, a:-H = 115, A/(a;) = (a;+ l)'-a;'= 1520875- 1481544 = 39331. This 
result is set on a line midway between the numbers corresponding to a;=114 
and a; = 115. In the same way the 2nd differences, A-, are obtained from the 
1st differences, the results being set on a line midway between the latter, 
and so on. 

It will be noticed that the 3rd differences of 3? are constant and that its 
4th differences are zero. It will be shown later that the mth differences of a 
function of the Jith degree are constant and that the differences of higher 
order are all zero. 

Example Z. — To determine the suooessive differences of sin (a a; + 6). 
Let /(a:) = sin (ok + 6). 



Then A/(a;) = sin (aa: + jo + 6) -sin (aa: + 6) 
=2 sin J a tocos (a a; + 6 + Jo to) 



= 2 sin iat«sin(oa; + 6 + 4oio + 7r) 



A^/(a;) = A{ 2 sin \aw ain (oa; + !> + Jow + ir)} 



=2 sin i o t« { sin (o a: + 1« + 6 + J o w + tt) — sin (a a; + 6 + J o io + tt)} 
= {2 sin iai«Poos(aa; + 6-f J2ato + ir) 



={2 8in4aioPsin (aa:-l-6 + J2a«) + 2 7r). 
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Similarly, it can be shown that 



^^/(x)={2 ami a wYain {a x + b + iSaw + Ww) 
and, more generally, that 

^"/(a;) = { 2 sin J a 2o}" sin [ax+b + ^naw + nn). 
This can be proved by induction in the following way :- 



A"+i/(a::) = A[{28in Jaw}" ain {ax +b + inaw + nv)'l 
={28in Jciuj}" [sin {ax + w + b + inaw + nir) 
-sin (ax + b + ^naw + mr)] 



={ 2 sin i a M)}''+i cos {ax + b + i(n + l)aw + mr) 



= {2 sin 4 a !«}"+' sin {a x + b + i {n+ I) a w + {n+ 1) v). 
But the theorem has been proved for 7i= 1, 2. Hence it is true always. 
Example S. — Find the first differences of the functions 

cos(aa; + 6); tan-' aa; ; a'"=. 
Example 4- — Find the »th differences of the functions 
— ; a' ; oos'(aa; + 6). 

X 

2. Properties of the Operator A. 

The operator A obeys certain of the fundamental laws of 
algebra, viz. — 

(i) The Law of Distribution. 

For A {/{x) + <^ (x)} = {/{x + w) + 4>{x + w)}- {/(x) + </. (x)} 
= {/{x + w)-f{x)} + {<i>{x + w)-<f>{x)} 
= Af(x) + Acl>{x). 
(ii) The Law of Commutation if the coefficients of the various 
terms are constants. 

For, if 01 is a constant, 

A [a .f{x)} = af{x + w) - af{x) 
= a{/{x + w)-f{x)} 
= a.A/{x). 
(iii) The Law of Indices. 

For A-- . A'/ (a;) 

= {A . A . A ... (r factors) . A , A . A ... (s factors)} /(a;) 

= {A. A. A ... {r + s factors)} /(a;) 

= A'-+'/(a;). 
The indices r and s are, of course, integers : the operator 
A" was defined in the previous section only for integral values of n. 
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Example 1. —If /[x) is a polynomial in x of the rath degree, then its nth 
difference is constant and its differences of higher order than ti are zero. 

Let f(x) = a„x!"+an-iX"-^+ ... +aiX + ao 

in which the coefficients o„, On—i, ... ai, oo are constants. 

Then A/(x)= {a„ {x + io)» + a"-' (x + !o)"-i + ... + oi (.r + w) + Oo } 
— {a„ a;" + a,i_i a;"-' + ... +aia; + ao} 
= 71 a„ X"-'' w + h„_-2 X"-- + ... +bix + ho 

where 6„_2, ..., ft], to are constants. 

We thus see that A/(a;) is a polynomial in z of degree (» - 1 ). Repeating 
the operation, we have 

AV W = {n a„ (a; + !«)'—izu + 6„_2 (a; + !«)»-=+ ... + 6i (a: -^ i«) + 60} 

- {n an x''-^ w + bn-2 x"-^ + ... + 61 a: + 60} 

=n(n-l)a„X"-^m' + Cn-3X"~^+ ... +CiX + Ca 

where c„—3 C], Cq are constants. 

Hence i?f(x) is a polynomial in x of degree (?i - 2). It is evident that 
by continued repetition of this process the degree of the function will be 
reduced to zero, and that we shall obtain 

^■^f(x)=n(n-\)(n-2) ... 3.2. \a,iio<< . 

i.e., the nth difference of the function is a constant. Obviously the next and 
all differences of higher order will be zero. This result is of great import- 
ance ; probably no other theorem will be applied as frequently as this. For 
instance, it enables us to complete a mathematical table in which there are 
gaps. This is illustrated in the next example. 

Example 3. — If 
/{a;) = 0-5563, /(a;f l) = 0-5682, /(a; + 3) = 0-5911, /(a;-l-5) = 0-6128, 

obtain /(a: + 2) and /(a; + 4) on the assumption that /(a;) is a polynomial of 
the third degree. 

It has been shown in § 1 that 

AV(a:) =/(a; + 2 w) -2f(x + w) +/{x). 
By repeating the process we get 

^*f{^)=f{3: + iic)-if(x + Sio) + 6f{x + 2w)-4/{x-\-iv)+f[x). 

But the given function is of the third degree. Hence for it we have 
A''/(a;)=0. Noting that w is unity we obtain 

/{a; + 4)-4/(a: + 3) + 6/(a: + 2)-4/(a:+l)-l-/{a:) = 
or, denoting /{x + n) by /„, 

/.-4/j + 6/,-4/i+/o=0. 
Similarly /s" */i + 6/3-4/, +/i=0. 
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Substituting tlie given values in these two equations, we get 
/4 + 6/2 =4 -0809 
/,+ /2=1-1819 
whence /2=0-5798 , /4=0-6021 , 

i.e., /(a; + 2) = 0-5798, /(x + 4) = 0-6021. 

Example 3. — If when a; = 0, 1, 3, 4, 6, the corresponding values of 
f{x) are 1, 1, 79, 253, 1291, find approximate values for f(x) when x=2 
and when ^,- = 5. Under what conditions will the result be accurate? 



3. Expansion oif{x-\-nw) for Integral Values of w. 

We shall now establish certain formulae which express 
J{x + nw) (where n is an integer) in terms of differences of 
various orders. These formulae will be shown in Chapter II. to 
be valid under certain conditions even when n is not an integer. 

It will be found convenient to introduce another operator E 
defined by the equation 

This new operator will obviously obey the same laws as A. 

Now EJ (x) = (1 + A) fix) =f{x) + A/{x) =f{x + w). 

Thus the effect of the operator E on the function f{x) is to increase 
its argument x by the given constant quantity w. Repeating the 
process we get 

E^J{x)=f{x + nw) 

and since the operator A may under the conditions of § 2 be treated 
as an algebraical quantity, this last equation may be written 

f{x + nw) = (\+Aff{x) 

=/(x) + ,AAy(a;) + „C2AV(a:)+ ... +AV(a;) (1) 

where „c,. denotes the coefficient of a' in the binomial expansion of 
(1+a)". The function /(a; + ww) is here expressed in terms of 
f{x) and its successive differences. 

It is evident that, since the effect of the operator E on f{x) is 
to increase x by w, the result of operating on f{x) with E~', where 
E~' is defined by the equation ££"^ = 1, will be to subtract w 
from X, 

i.e., E-^f{x)=f{x-w). 
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Hence the above formula (1) may be transformed as follows : 

./ {X + nw) =/(x) + „c, ^f{x) + (1 + A) U-' {„c„ A'fix) + „c, A'f{x) 
+ ... +A"/W} 

=/(x) + ,A Ay (a:) + (1 + A) {„c., A= /(a; - w) 
+ „c,AY(a;-M;)+ ... +A"f(x-w)} 

=/ (a;) + ,.c, A/(a;) + „c„ A=/(a: - w) + „+,c, A'/{x - w) 

=/(a;) + ,A A/(jj) + „C2 A-f{x ~w) + „+iC3 AV'(a! - w) 
+ {\+A)E-' {„+,c. Ay (x - M)) + . . . } 

=f{x) + nCi A/(a!) + ,fi„_ A'f{x -w) + „+iC3 A'/(a; - m;) 
+ (l+A){,.+AAV(a:-2w)+ ... } 

=f{x) + ,fi, A fix) + A A' fix -w) + „^,c. A' fix - w) 
+ u+,c. A-' fix - 2w) + „^,c. A' fix - 2m;) + ... 

Operating next with (1 + A) £!~^ on the 7th and subsequent 
terms, then on those beginning with the 9th term, and so on, we 
obtain 

fix + nw) =/(x) + „Ci Ay (a;) + „Co A=/ (a; -w) + ,,+iCs AV(a! - w) 

+ n+ic, AV(x - 2«;) + „+,c, A' fix -2w)+ ... (2) 

the formula ending at the term A-"-'/(a: - n - I w). 

In the same way it may be shown that 

fix + n w) =fix) + ,fij A fix -w) + „+iC„ Ay (£c - w) 

+ n+A Ay (a; - 2w) + „+,c, A^ fix - 2w) 
+ „+x,A'fix-3w) + ... (3) 

the formula ending at the term A-"fix - n w). 

Changing x into x + w and n into ?i - 1 in (3) we have 
fix + n w) =/(x + w)+ „_,ci A fix) + „C2 Ay (x) + ^c. Ay (a; - w) 

+ „+iC4 Ay (k - w) + „+iC5 A' fix - 2m;) + . . . (4) 
the formula ending at the term A-"fix - n- 2 w). 
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Adding (2) and (3) we obtain 

J.- 1 \ ri s ^f{^) + ^f{x-w) n ^, _ 

J{x + n w) =/ (a;) + ,fi, -^ ' ' ^'^ ^ '- + — ,fi, A7 (x - w) 

^^f{x-w) + l\\f(x-1w) 



2 



+ —n^,c,^*f{x--2w)+... (5). 



n 
T 



Similarly, addition of (2) and (4) gives 

Ax^n^)j:^^±^^^^^Hn-k)^m 

A'f{x) + A'/{x-w) n-i ^,^, 
+ „c, -"^ ' ^-^^ ^ + —^ „c„ Ay (as - w) + . . . (G). 

Formulae (1), (2), (3), (4), (5), (6) may be compared with the expansions 
introduced in Chapter II. under the names Newton's, Stirling's, etc.. 
Formulae of Interpolation. 

Again, from the definition of the operator U, we see that 

A"/(x) = (^-l)"/(a;). 
Expanding the operator in this last relation we obtain 
A"/{x) = E"/{x) - „c, E"~'f{x) + ,fi, m'-'f{x) ... +{- l)"/(,'«) 



=f{x + nw) - ,fii/{x + n-lw)+ ,finj\x + « - 2 w) 
+ ... +(-l)"/(a;). 

In this result we have the ?tth difference of the function f(x) 
expressed explicitly in terms of f(x) and its successive values. 

4. Relations between the Differences and the Dif- 
ferential Coefficients of a Function. 

The second example of § 1 may have reminded the reader of a 
somewhat similar result in the differential calculus, viz., 
J] 

-r- sin {ax-\-h) = a" sin (ax + b + imr). 
ax 

This result may be derived from the example in question. For if 
we put w = Ax, where Ax denotes a small increment in x, and. 
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after dividing both sides of the equation by (A a:)", make A x tend 
to zero, we get 

d"f(x) J /-2 sin A a . A a;-) " . , , , , . , ^^ 

dx" l-^\ ^x ) ' ^ ^ 

Ax->Q 

= a" sin (a X + 6 + J Ti tt). 

This naturally suggests that there may be some relations connecting 
the successive differences of f(x) and the successive differential 
coetiicients of the same function. That such a relation does exist 
we now proceed to show. 

Applying Taylor's Theorem* to the function /(a; + w) we obtain 

t(r w'" 

f{x + w)=f{x) + wf' {x)+ -f"{x)+ ... + ^/""'(«')+- 

or, since f{x + w) -fix) = A/ (a;), 

hf{x) = wj'{x)+'^^r(x)+ ... +^/M(aj)+... 

which is a relation of the kind sought for. 
Again, by Taylors Theorem, 

f(x + xo) ^J\x) + wf {x) + -f"{^x)+ ... + — /"»' {x)+... 

f{x + 1w) =/H + Iwf (x) + L_L/" (^) + . . + i__L /(".I (a,) + . . . 

and generally 

J\x + n w) =f {x) + n wf (x) + L_l /" {x) + ...+ ^— f- /•"" {x) + ... 

Hence 

A"y(x) =/{x + n w) - „c, f{x + n-\ w) + „c^'(x + n-2w) 

- ... +{-irf{^) 

w'" 
= 2 — {»i--„Cj(w-l)'" + „c,(n-2)'"- ... }/i'»'(x). 



* It would be inopportune at this stage to introduce a discussion of the 
questions of convergence which arise in connection with the application of 
Taylor's Theorem. We shall assume a convergence sufBcient for our 
purpose. 
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But 

n"'-,fii{n-iy" + „c„(n-2)"'+ ... 

= if m-Cn* 
= n ! ii T)i = n 

= — ?i(3n+l)(n + 2) ! ifm = ?i + 2 
etc. 
Making use of these results in the expi-ession for A"/' (sc) we find that 

15 n-(n + 1) ,..,„,,. 

Corollary 1. — Since the wth differential coefficient of a function 
f{x) of the Tith degree is a constant, and the derivatives of higher 
order are zero, we see again from this expression for i\^'f{x) that 
the nth difference of a polynomial of degree n is a constant, and the 
differences of higher order zero. 

Corollary 2. — By reversing the series for ^"/{x), we obtain an 
expansion for the nth differential coefficient of /(x) in terms of its 
differences, viz. 

w"/'"^ {X) = A" fix) - ^ A»+y (X-) + 'lil^^ A"+V(»;) 

_ 15n(...2)(. + 3 )^.,.^^^^^_ 

Corollary 3. — If the wth differences of a function are constant, 
then the function is a polynomial of the nth degree. 

For since the nth differences are constant, the higher differences 
are zero, and therefore the nth differential coefficient is constant, 

* Chrystal's Algebra, Part ii., Chap, xxvii., § 9. 
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and the differential eoeflScients of higher order zero. Thus the 
function, when expanded by Taylor's Theorem, takes the form 

f(x) =/(^) + (X - X)f (X) + ^-^^V" (^) + • • ■ 

and this is of the wth degree in x. 

Corollary 4. — Since e*=l+T— + 7r-;+ ••• -I ; + ••• 

\\ Z\ n\ 

we see that the series which represents f{x + w), viz. 
f{x) + wf{x)+-J"{x) + ... 

may be symbolically represented as 

e'^fix). 
a 

Hence Ef{x) = e "" .f(x). 

MISCELLANEOUS EXAMPLES. 
(Assume in the following exercises that lo is unity, unless otherwise indicated ). 

1. Find the first diflferenoes of the functions 

tan a x ; log x; (a"^ - a")! (a-\). 

2. Prove that 

sin 8 



A tan xS= ■ 



' oosafl oos(a;+l)fl ' 

3. Prove that 

1 1 e 

2' tan— 

4. Find the ?i"' differences of the functions 

x + l 
— — — — - . j;>n+n . sin ax cos bx. 
x^ + 6x 

5. Find the reth difference of each of the factorials x(x + a)(x + 2a)... and 
, there being m factors in each factorial and the increment 



x(x + a){x + 2a)... 
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of X being a : and show how to develop any polynomial in a: in a seiies of 
factorials. 

In particular show that 

ar" =a:+7a;(a;-l) + 6a;(a;-l)(a:-2) + K(a;-l)(a;-2)(a;-3). 

6. If /(a:) = a3* + 62"= prove that 

AV(a:)-3A/{a:) + 2/(a:) = 0. 

7. Prove that 

Ak" - 26?x" + 3 A^rE" - =(a;-l)" - (a; -2)". 

8. Establish the formula for tt e 7ith difiference of the product of two 
functions 

A''{/l.^■). 0(.-c)} = A"/(.^•)•0(a.■) + m A"-i/(.^+ 1) • A 0(a;) 

, n(n- 1) 
+ -^yi- '^"-•'/(a; + 2) . AVW + . . . 

9. Show that 

<t> {E) a«/(a;) = a' (a E)f(x) 
and find the value of 

AV(.r) + — (a - COS fc) A/(.r) + (" 1 - -?- oos 6 + 4r)/ W 
a \ a a- / 

cos t a; 
where /(a:) =— ^j— . 

10. The record of exports for the year 1813 was destroyed by fire. Make 
an estimate of their value given that the values of the exports for the years 
1810, 1811, 1812, 1814, 1815, 1816 were 48, 33, 42, 45, 52, and 42 million 
pounds respectively. 

11. It is asserted that a quantity, which varies from day to (?ay, is a 
rational and integral function of the day of the month, of less than the fifth 
degree, and that its values on the first seven days of the month are 30, 30, 
28, 25, 22, 20, 20. Examine whether these assertions are consistent. If so, 
assume them to be true, and find (1) the degree of the function, (2) its value 
on the sixteenth day of the month. 



[CH. II 



CHAPTER II. 

FORMULAE OF INTERPOLATION. 

5. Introduction. 

The problem of determining the value of a function for 
some intermediate value of the independent variable or argu- 
ment, when the value of the function is given for a set of discrete 
values of the argument, is called Interpolation. When the function 
obeys some definite law the approximation to the true value may be 
made as accurate as may be desired : but if, as often happens, 
the rigorous analytical formulae which occur are very com- 
plicated, they will be of little value from the point of view of 
computation. It will then be necessary to use a method analogous 
to that which has to be employed when the properties of the 
function are entirely unknown, save for the values (usually obtained 
as a result of experiment or observation) which it has for a certain 
limited number of values of the argument. It is to this case that 
our attention will be directed. 

When the observations are given for values of the argument 
equally distant from each other, and are in arithmetical progression 
so that they can be plotted on a straight line, the possibility of 
interpolating is obvious. But, in general, the observed values 
do not follow any apparent law, and unless some assumptions 
are made regarding the form of the function, the problem is 
indeterminate. It is, however, customary to obtain the entry for 
values of the argument not far distant from each other ; and when 
this is done, in most practical problems it is found that all the 
differences at some particular order become sensibly zero. In such 
cases the observations may be represented with considerable 
accuracy by a polynomial in x of the form 

2/ = a,, + flSi a; + 02 a? + ... + o[„ a;". 
This is the equation of a parabola of the wth order, on which ' 
account the method now to be described is sometimes called 
Parabolic Interpolation. 
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In some cases it is advisable to calculate the function not 
directly, but in combination with some other function. The 

calculation of the values of the probability integral e-i*" eft 

furnishes an example of this kind, since, when the argument becomes 
large, the integral and its differences become inconveniently small 
and irregular. In this case it will be found convenient to compute 

the function e^''' I e-i^'dt. At other times some function of the 



■j>-"' 



given function is calculated ; for example, if the observations 
formed a geometrical progression, then we should prefer to make 
the interpolation on their logarithms rather than on the quantities 
themselves. 

Wallis (1616-1703) may be said to have laid down the principle of inler- 
polation in his approximation to the area of the circle 2/ = (.t- .r^)'/^. In this 

work he first assumed that the value of the integral I {x - x^)^ ' dx might be 

Jo 
taken as the geometric mean of the values of • 

I {x - x'^)" dx and { {x-x^)dx, i.e. of 1 and J. 

Subsequently he saw that this was not exactly true, and that the value of 

I {x- x'^yl- dx must obey the law expressed by the series 
Jo 

[ (x- x^f dx, \ {x-x')dx,..., I {X- x^)« dx ; 

Jo Jo Jo 

this was equivalent to interpolating the value of the integral under con- 
sideration. 

6. Lagrange's Formula of Interpolation. 

We shall first show how to find a polynomial of degree (n - 1) 
at most, which has the same values as some arbitrary function /(x) 
for n distinct values of x, say a; = ai, a^, ... a„. This polynomial 
will be represented by the parabola of order (n-l) already men- 
tioned, and may be expressed in the form 

/(as) = Cj + Ci as + C2 a;- + . . . + c„_i x"-'. 
But it will be found more convenient to write it, as is always 
possible, in the form 

/(«)= 2 ^(«-«i) ■■• (as - a,_i) (as - a,.+i) ... (x-a„) (1) 
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When X is put equal to a,, all the terms on the right-hand side of 
(1) vanish except the term containing A,: at the same time y (a;) 
becomes f{a^. Hence 

^r =y (O / («r - ail) . . . {O'r ' «.-l) (^r " ^r+l) • ■ • {^r " «„) 

and equation (1) becomes 

/(a,) = 2 (a^ - «i) • ■ ■ (a^ - «.-i) (a^ - «.+i) ■ • ■ (a? - a„) ^^^^^ ,3^ 

r = l (»r - »l) ■ • • (»r - «r-l) (»r " "r+I ) ■ • • ("v " ««) 

This formula expresses f{x) for all values of x in terms of its 
values at the n points a,, o,, ... a„. If f{x) is a polynomial of 
degree less than n, the formula gives its accurate expression. But 
if y(a;) is not a polynomial of degree less than n, then it is no longer 
possible (in default of other information) to determine f{x) from a 
knowledge of its values at the points a^, a^, ... a„; for if any one 
such function were known we could obtain another by adding to it 
any function which is zero at the points a^, a,, ... a,„ e.g., we could 
add the function (a; - a,) (x - a,) ... {x - a„) 4>{x) where ^(a;) does 
not become infinite when x = a,, a^, ... a„. If, however, f{x) is a 
function of which nothing is known except that it takes certain 
definite values corresponding to certain definite values a^, cio, ... a„ 
of the variable, then the simplest assumption which we can make 
regarding its form is, that it is the function of degree (^ - 1) which 
satisfies the given conditions ; and this function is given by 
equation (2). The expression on the right-hand side of (2) may 
therefore be taken as the expression for f{x) for all values of x in 
the interval considered. It may be written in the form 

fix) = ^ -, - ,, , , /(a,) 

^=1 (a; - a,) <^ (a,) 

where <^ [x) = (a; - a{) (x - a^) ... (x- a„). 

This is known as Lagrange's Formula of Interpolation. The 
advantage which it has over other formulae for the same purpose 
is that it is in a form suitable for logarithmic computation. 

This latter form of Lagrange's formula may be obtained directly from the 

formula previously given. It can also be obtained very readily by the 

Method of Partial Fractions. 

fix) 
For consider the proper fractional function ^ , , where 

(t>{x) 

<t>(x) = (x-a{) {x-a^... (x-a„). 
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Decomposing it into partial fractions we have 



+ ... + ——-. 
^{x) X -a^ X- a^ X - o„ 

But, by the general theory of partial fractions, the coefficient Ar is given by 

filr) 



Ar = 



<t>'(a-r) 



Hence 




















f(x) 
<t,(x) 




/(»i) 


)^- 


{X- 




■■■^W- 


/( a.. ) 




(X- 


-<h)<t>'(<h. 


- a„ ) ip'{a„ ) 


which 


gives 




/{=) = 


n 

2 




fix) 


f{ar). 






(X 


-ar)<t>'{ar) 





Example 1. — Construct a polynomial of the third degree in x which 
shall have the values -5, 1, 1, 7 when the values of x are -2, - 1, 1, 2 
respectively, and determine its value when x is 9. 

Let f(x) be the required function. Substituting the given values in 
Lagrange's Interpolation Formula we have 

(.-g + 1) (a;-l)(a:-2) ,- . _, ... _, ,,. _. 

f(x) — j_2+i)(_2-i)(-. 2-2)^"^'"'' '-ij-fliw-i-Hf-i -9^^ '"^) 



(a; + 2) (x 


-l)(x- 


■2) 




(-l + 2)(- 
(x + 2) (x 


1-1)(- 

+ 1) (x 


-1- 


-2) 
1) 



(X + 9) {X + 1) (a: - 2) , , ^^ , \-^ -r ^; ^u, -r x, ^.^ - ., . _ 

(I + 2) (1 + 1) (1 - 2) ^ (2 + 2) (2 + 1) (2 - I) 

= i%(a;+l)(a;-l)(a;-2) + J(a; + 2)(a;-l)(a:-2) 

-Ha: + 2)(a; + l)(a;-2) +/i (.-k + Q) (a;+ 1) (a:- 1) 

= a^-a; + l. 

From this it follows that /(^) = 25. But in practice it is not necessary to 
obtain the algebraic form of f(x) if we merely require its value for some 
value of X. Direct substitution in the formula will give the result much 
more quickly. In the question under consideration 

(f+l)(#-l) (^-2) (t + 2) (j-1) (t-2) 

.A(|)-(_2+l)(_2_l)(-2-2)>~°'+ (-1 + 2) (-1-1) (-1-2) 

(t + 2) (§+1) (t-2) (4 + 2) (j+1) (4-1) 
+ (1+2) (1 + 1) (1-2) + (2 + 2) (2+1) (2+1) ^ + ^' 

_ 26 7i 35 1245 

Example 2. — Find the minimum value of the polynomial of the second 
degree which has the following values 

X 0, 1 , 5 
f(x) 11, 5, 21. 

a. 2 
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Substituting these values in Lagrange's Interpolation Formula we have 

(a;-l)(a;-5) x(x-5) x(x-l) 

/W = -1.-5 -It + TT^-S + -574^- 21 

.-. /'(:») =-V- (2 a; -6) - |- (2a:-5) + |^ (2a:- 1). 

For a minimum turning value /' {x)=0 

i.e., 44(2.-i;-6)-25(2a;-5) + 21(2a;-l) = 

80 a; -160 = 
x = 2 

Hence the minimum value oif(x) corresponds to the value 2 of the argument. 
The corresponding value of f(x) is 3. 

Examples. — If /(x) is a function of x which assumes the values 20, 35, 56, 84, 
when X has the values 0, 1, 2, 3 respectively, express f(x) in terms of x, on 
the supposition that A*/(a;) = 0. Also find an approximate value of x if 
/(a;) =25, third differences being neglected. 

Example 4- — If four consecutive values of a function corresponding to the 
values 1, 2, 3, 4 of the variable be 601, 777, 902, 999 respectively, find the 
value of the function when the argument has the value 1 '5. 



7. Periodic Interpolation. 

When there is reason to suppose that the function is periodic, it is better 
to assume for /(a;) an expansion in terms of periodic functions, say, 

/(a;) = Oq + Oi cos x + ttg cos 2 a; + ... + a„ cos n x 
+ bi sin ai + ftgsin 2a;+ ... +6„sinna;, 

the number of unknown coefficients corresponding to the number of observa- 
tions given. When the observations are given for values of the argument 
spaced at equal intervals, the problem becomes one in Fourier's analysis and 
will be found treated in the tract on that subject. But when the values of x 
are not in arithmetical progression we may use a formula which is analogous 
to Lagrange's formula of interpolation, and which Gauss has shown to be 
equivalent to the above expansion for /(a;). 

If /(x) is given for a;=a;Q, .a;,, aij, ... x^n the corresponding formula is 

fi~\ -^ sin ^ (a: - a:;,) ... sin |(a;-a ;r_i)sin 4 (a;-a;,+i) ... sini(a:-a;2n) 

J\X) — ^ J\Xr]' 

r^o sin J {xr - ail) ... sin 4 (x,. - a;r_i) sin J {x,- - Xr+i) ... sin J (Xr - X2„) 

8. Notation. 

Hitherto the points a^, a^, ... a„ have been taken to be any 
points whatever. But in the most important class of cases it is 
possible to obtain the entry for values of the argument equally 
distant from each other. In what follows we shall assume that the 
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arguments are spaced at equal intervals, and shall use the notation 
of the following scheme which for many purposes * is more con- 
venient than that used in Chapter I. 

A8 



Argument. 


Entry. 


A 


A2 


A3 


A< 


A» 


a-3w 


/-3 


«/-? 










a -'2w 


/-2 


Z^S/.a 


sy.; 


2 








M/-? 


S/-0 


/^sy-; 


; sy.; 






a - w 


/-I 


p-¥-i 


sy-1 


l"Sy-; 


1 sy.i 






F/-i 


S/-i 


f^sy-i 


sy.. 


\ f^^V-i 


sy-i 


a 


fo 


/^S/o 


sy 


/'Sy 


sy 


A'sy 




M 


Vj 


/^sy 


syj 


/^syj 


sy 


a + 7.v 


/i 


y^S/ 


sy 


/^sy 


sy 






i^fi 


S/i 


Msy 


sy 






a + '2w 






sy 








a + 3w 


A 












where 


V-i 


=/»-/- 


1 


/^/i = i 


(/o+/i) 






8A- 


= /i-/o 


J 


"S/„=| 


(S/i + 8/. 


i) 






etc. 




etc. 





sy 



9. Newton's Formula of Interpolation. 

If we assume that the diiferences aj - a^, a^- a^, ... in Lagrange's 
interpolation formula are all equal to w, and that Oj in the formula 
corresponds to a in the scheme, then the quantities /(a^) take the 
following forms : — 

/K)=/.=/i + s./« =/o+2S/j + sy 
f(a...}-MrSA ^ ^-V^ sy -H !±z^(!:z^ sy, ^ ... 

* Especially in connexion with central-difference formulae. 
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Lagrange's Interpolation Formula then becomes, if the number 
of given points is denoted by k, 

^ (x - gg) (x - a,) (x-at) ...{x- a,) 
(ai - a,) (tti - as) {a^ -a^) ... {a^ - a„) 

{x-a^){x-as){x~ai)...{x-a^ .. 

(«2 - «i) («2 - «3) («2 - (ij) • • • («2 - «l) 

(«3 - «l) (*3 - ftj) («3 -di) ■■■ (^3 - %) 
+ 

(x - a,) ... (a; - a,) (a; - a^+o) ■■■ {«>- «*■) 



+ 



(«,+l - Ol) ■■• («,.+! - "r) («,■+! - «^+2) ••• («--+l - «s) 



+ 

where ai = a, a2 = a + w, ..., a^ — a + {k-l)w. 

This is an expression of degree (A; - 1 ) in a;. In it the coefficient 
of /o is unity when x = a^, tiai • • • '"■a ^^'^ hence is unity always. 
Again, the coefficient of S/j is at a^, 1 at a^, 2 at a^, ... r at 
a,+i, . . . , and as this coefficient must be a polynomial in x of degree 

less than k, it must be no other than . In fact the terms in 

w 

8/ in the above series simply constitute the expansion of 

ryt /Tf on yji 

8/. by Lagrange's formula. We shall denote by n, so 

that M is a real quantity, not necessarily an integer. The 
coefficient of 8/. is thus n. Similarly the coefficient of S^/^ is 

r (r — \) 
at O], at a,, 1 at Oj, ... — - — — at a,^, . . . , and in the same way 

we see that this coefficient is — ^-r-j — - for any value of n. The 

other coefficients can be determined in the same way, and we 
thus finally have Lagrange's formula transformed into the equation 

/•/ , \ J- . s:r , '^('*-l) S2J- ^(»i-l) ('^-2) ., . 
f{a + nw)=fo + nSfi + ^ , ^A + jy -' Sy^ + ... 
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This is called Newton's Formula of Interpolation* In it the 
differences are those which occur along a line parallel to the upper 
side of the triangle in the scheme. It may on this account be 
termed a Forward Difference formula of interpolation, and will be 
found very useful when the value of the function near the beginning 
of a limited series of observations is required. 

It must be understood that the function furnished by Newton's 
formula is the polynomial of lowest degree which has the given 
values at the k given points a, a + w, a+2w, ... , a + (k-l)w. 
and the series on the right is a terminating series. If the given 
values at the k given points are the values of some given function 
<f) (x) at the points, we cannot conclude that 4> (x) is represented by 
Newton's formula unless <^ (x) is the polynomial of lowest degree 
which passes through the points. 

Suppose, now, that the values of a function are given not 
merely for a finite set of values of its argument, but for the infinite 
set 

x = a + r w (r = 0, I, 2, 3, ... ad inf.). 

Then we can, by Newton's formula, construct a polynomial of 
degree (k-l) which has the given values for the first k of these 
values of the argument, where k is any positive integer : and it 
may happen that, as k increases indefinitely, the series on the right 
in Newton's formula becomes (for some range of values of a;) a con- 
vergent series, having for its sum a function /(x). The function 
J'{x) will now have the given values for each of the values a + rw 
(r = 0, 1, 2, 3, ... ad inf.) of the argument. There are, however, an 
infinite number of functions which satisfy this condition and f{x) 
is that one of them which is the limit, when ^ -> oo , of the polynomial 
of lowest degree which has k of the given values. 

An expression for the difference between the value of the 
function and the sum of the first {r + l) terms of the limiting 
form of this polynomial will now be obtained. Let 

nin-l) ... (n-r + l) . , 
f{a + nw)=f + nSfi + ... + -^ '—^ •' S7_. + R. 

Then 

m(m-l) ... (m-r+l) 
f{a + mw)-f-mS/i - ... - -^ '-^ ' B'f^ -R 



Cf. § 3, Formula (1). 
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will be zero for m = 0, 1, 2, ... r and n : its derivate of order (>•+ 1) 
will therefore be aero for some value n' of »n lying between the 
least and the greatest of the quantities 0, r, n, and hence the 
residue is given by 

R = „C,+, m;''+^ /"■+" (a + n' w). 

It may be remarked that Newton discovered the binomial coefficients 
when studying the expansion in aeries of a function by interpolations. He 
first of all considered the expansions of the expressions (1 -a;^)", (1-a;-)^ 
(1 -aj'p, ... and deduced from them the expansions of (1— a;^)i, (1 -a;^)^, ... . 
Then, by analogy, he obtained the expression for the general term in the 
expansion of a binomial, and thus established the binomial theorem. 

Example, 1. — From the following table of values of the function f(x) 
determine the value of /(3'5). 

X f(x) A A2 A3 

1 0-208460 

29242 

2 0237702 - 213 

29029 - 27 

3 0-266731 -240 

28789 -26 

4 0-295520 -266 

28523 

5 0-324043 

In this example a=\, n=^, w=l 

/,=/{!) = 208460. 
Hence /(3-5) =/„ + ^ 5/^ J- sy, + A' S^f. . 

=0-208460 
14621 
27-2 
=0-223106. 

Example S. — Corresponding values of x and /(a;) are given in the following 
table : — 

X f(x) 

4 54-5982 

4-1 60-3403 

4-2 66-6863 

4-3 73-6998 

4-4 81-4509 

4-5 90-0171 
Obtain the value of /(4-05). 
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A certain amount of discretion is necessary in the application 
of Newton's formula. Let us assume, for instance, that the 
given observations are finite — say seven — in number, and 
correspond to the values a- 3w, a-2w, ... a + 3iv of the 
argument as in the scheme, §8. Ji a + nw lies between a-3w 
and a-'2w, then Newton's formula iorf{a + nw) will involve the 
quantities /_3, S/_j, 8y_,, SV_^, 8y_„ S'f_^ and 8y„, all of which 
occur in the scheme. Thus a good approximation to the value of 
f{a + nw) may be obtained by means of the formula. But if 
a + nw lies between a + 2w and a + 3w, then in the Newbon's formula 
for /{a + nw) the only known terms are the first two depending 
on /2 and S/s. In general, these are quite insufficient to determine 
the value of /{a + nw) to the necessary degree of accuracy. It is 
therefore expedient in this case to find a formula involving a 
larger number of known quantities. A suitable one may be 
obtained by replacing w in the original formula by -w. Then, 
since 

we must replace S/l by 

fo-f-i or -¥-h- 
Similarly, since 

we must replace B^J\ by 

In the same way we can determine the functions which must 
take the places of SySi, S*/^ Newton's formula then becomes 



./■(« ■ 






2 ! • ' 3 ! 



This expansion involves differences along a line parallel to the 
lower side in the scheme, and it may therefore be regarded as 
a formula for Backward Interpolation. When a + nw lies between 
a + lw and a + 3w, the formula will depend on the quantities 
A, 5/s> Syo, SVs, Sy^, S'/i, and 8y„, all of which are known: it 
is therefore very convenient when the value of fix) is required for 
values of x near the end of the series of observations. 
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Example 3. — From the following table of values of /(a:) determine the 
values of /( 104 -25). 



X 


f(x) 




A 


6? 


A' 


101 


1030198 




30906 






102 


1061104 




31518 


612 


6 


103 


1092622 




32136 


618 


6 


104 


1124758 




32760 


624 




105 


1157518 










In this example a 


= 105, n = 


0-75, 


w= 1 








/o=/(105) 


= 1157518 






/(104-25) 














0-75 (0-75 


-•)^ 


0-75 (0-75 


-l)(0-75- 


-2) 


=/o-0-15Sf_ 


^+ ., 


s^y-i- 


3 ! 





-3+... 

=/o-0-755/'_j -0-093756-/_j -0'039062o5Vs + ... 
= 1 157518 - 0-75 (32760) - 0-09375 (624) - 0-0390625 (6) 
= 1157518-24570 
58-5 

0-234375 
= 1132889-265625. 
Example 4. — Using the data of Example 2, find the value of /(4-45). 

10. Stirling's Formula of Interpolation. 

As Newton's foi-mula involves differences which lie along a 
slanting line it will frequently not be suitable for calculating the 
value of the function for a value of the argument near the middle 
of a series of observations. By a suitable transformation of this 
formula, however, we may derive another involving differences 
which lie on a horizontal line ; such a formula is termed a Central- 
Difference Formula. 

Since 



we have 
Hence 



8V; = 8Vi+/^8y„+i8v„. 
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In the same way the higher differences which occur in Newton's 
formula may be determined in terms of the differences which occur 
on the same horizontal line as /J,. Substituting these in Newton's 
formula we obtain 

=y„ + « M 'Vo + ^ 87„ + ^g, // 87„ + ^^, 8Vo 

M(n--l)(/t^-4) 

+ gl /^ 87„ + ... 

This is known as Stirling's Formula of Interpolation. It should 
not be used for determining the value of a function for values of the 
argument near the beginning or the end of a series of observations, 
as it would not then be possible to obtain a sufficiently high order 
of differences. 

Example 7.— Using the data of §9, Example 1, determine the value of 
/(3-5). 



X 




/(^■) 




A 


A' 


A3 


1 




0-208460 




29242 






2 




0-237702 




29029 


-213 


-27 


3 




0-266731 




(28909) 
28789 


-240 


(-27) 
-26 


4 




0-295520 




(28656) 
28523 


-266 


(-26) 
-26 


5 




0-324043 




28-231 


-292 




6 




0-352274 










this case 


u 


= 3, w= 1, 


n = 


--h 








/o 


= /(3) = 0-266731, 








/(3-5) 


= /o + 4mS/o 


+ ^ 


^Vo-A5V„ 










= 0-266731 














14455 














2- 


30 












= 0-281158. 











Example ^.-Using the data of § 9, Example 2, find the value of /(4-21). 
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If n is greater than 0'5, it is better to use backward differences. 
The corresponding formula may be obtained by merely changing 
the sign of w in Stirling's Forward Difference Formula, viz. : — 

f{a - n IV) =y; - n /. 8/, + ^ S% - "^^^ M S^o 

n-(n=-l) n (71^ - I) (n' - i) 
+ 4T~ Vo 5"] — — ^5y„+... 

Of course, either formula may be applied whether n is greater or 
less than 0-5 : thus, when a check is necessary, it will be found 
better to computey(a;) by the formula not already used. 

Example 3. — Use this Backward Difference Formula to check Example 1. 
Here a = 4, i« = 1, n = ^ 

/„=/(4) = 0-295520 
/(3-5) = /„ - 4 M S/o + i S=/o + tV «Vo 
= 0-295520 - 14328 
33 
2 
= 0-281157. 

The agreement with the former result is perfect when we take into 
consideration the fact that the last figure is forced in both cases. 

Example 4. — Using the data of § 9, Example 2, find the value of/(4-29). 

11. Gauss' and Bessel's Formulae of Interpolation. 

Another central difference formula may be obtained by trans- 
forming Stirling's formula. 

Since /^8/u = Sy,--iS7; 

fiS'/o = SY^-^S*Jl etc., 

we see that Stirling's formula may be written in the form 

, (n + -2)(n+l)n(n-l)(7t -2)„, , 
+ gT 84+ - 

In this result, which is frequently known as Gauss' Formula, 
the /u-terms do not occur. The differences which occur in this 



5v;=i"8=y;-i8V' etc. 
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formula are those which lie on two adjacent horizontal lines of 
the scheme in § 8. 

Again, since /; +/; = 2 /j./^ 

we have 

and hence 

(n + 2)(7. + l)n(r.-l)( 7i - 2) 

+ 5I s/^+ ... 



3! 
(»i+1)m(w- l)(»i-2) 



4! 



/^s% 



^ (n + l)n(rt-l)(n-2)(n-^) g^^. ^ 

This is frequently known as Bessel's Formula. In it the differ- 
ences are those which occur on a horizontal line mid-way between 
the lines on which yj and /J stand. 

The formula most suitable for the construction of mathematical 
tables, when differences of higher order than the fourth are 
neglected, is obtained from the first of these on putting 

When this substitution is made we get 
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Bessel's formula will be found particularly suitable when n is 
equal to ^, for in this case the formula becomes 

in which only even differences occur. 

As with Stirling's formula, none of these results should be 
employed when the value of the function is required for values 
of the argument near the beginning or the end of a series of 
observations. Stirling's formula will be found more convenient 
when n is very small, or when the number of observations given 
is odd : Bessel's when n is nearly equal to a half, or the number 
of observations even. 

The backward difference formula corresponding to Bessel's 
forward difference formula may be obtained in the same way as in 
Newton's and Stirling's cases. But, except as regards sign, the 
coefficients in the two formulae will be the same, and therefore the 
one cannot be used as a check on the other. When a check is 
required to a result computed from Bessel's formula, it is better 
to compute it again using Stirling's formula. 

Example i. —The values of a function for the values -2, —1, 0, 1, 2, of 
the argument are respectively 0-12569, 0-17882, 0-23004, 027974, 0-3-2823. 
Using Gauss' formula, show that when the argument has the value 0-4, the 
value of the function is 0-25008. 



X 


/{^) 


A 


A'^ 


A3 


-2 


0-12569 


5313 






-1 


0-17882 


5122 


-191 


39 





0-23004 


4970 


-152 


31 


1 


0-27974 


4849 


-121 





2 0-32823 

In this example a = 0, 9i = 0-4, to = 1 
/„=/(0) = 0-23004. 
Hence /(0-4) = /„ + 4 S/^ - 0-12 S'f^ 0-056 S^f^ 

= 0-23004 
1988 
18-2 

= 0-25008. 
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Example 2. — Using the data of § 9, Example 2, find the values of /(4-21) 
and /(4 -29). 

12. Evaluation of the Derivatives of the Function. 

When, as often happens, the values of the differential co- 
efficients of a function which is given by a series of observations 
are required, use may be made of the synjbolic equality established 
in § 4, Cor. 4. 

For, since 

e'"^/(x) = Ef{x) 
we have 

w^.f{x) = ]ogE.f(x). 

CiOC 

The method will be evident from the following example : — 

Example 1. — The values of a function corresponding to the values I'Ol, 
1-02, 1-03, 1-04, 105 of the argument are 1-030301, 1-061208, 1 092727, 
1-124864, 1*157625 respectively. Find the value of the first differential 
coefficient of the function when the argument has the value 1 -03. 

Here w = 0-01. Let /„ = /(l-Ol). 

Then w j^f( 1 -03) = log £? . / ( 1 -03) 

= log^./2 

= \ogE.lPf, 

= log(l + A).(l + A)2./„ 

/ A^ A3 \ 

= (^- T+ T -■■ )-(l+2A + A'-')./„. 

Neglecting dififerenoes of higher order than the 4th, this gives 
dx- 



0-01 -4, /(I -03) = (A + ? AH 4 A^ - A AVo 



= (^-lE+%E^-^&)f, 

= l(/3-/l)-tV(/4-/«) 

= 0-042437-0-010610 
= 0-0318-27 
/'(l-03) = 3-1827. 



Example «.— The values of a function corresponding to the values 1-5, 2, 
2-5, 3, 3-5, 4, of the argument are 4-625, 2-000, -0125, -1, 0-125, 4000, 
respectively. Find the value of the first differential coefficient of the function 
at the points where the argument has the values 2-5 and 3-0. 
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The method is perfectly general. The same result may be 
obtained by differentiating any of the formulae of interpolation 
already established. For example, if we differentiate with respect 
to n the formula of Stirling, 

/(a + n m;) =/„ + «/. S/„ + Ij- Sy„ + !il!^) ^ 63/„ 

we have 

wf (a +nw) = iJ.S/o + nS^ f„ + '}{3n- - l)iJ.S''f„ 
+ 3-V(2»i-'-n)5y„+ ... 
w'f" (a+nw)^ 8y„ + n /. 8V„ + yV (6 rr - 1) S^o + ■ ■ • 

Putting n = these become 

w/' (a) = /i 8/o - 1 /i Syo + sV /^ S'/o - . . . 

wV" (a) = svo - tV syo + ^ syo - . ■ • 



These formulae are of use in determining the turning values of 
a function, the values of the argument corresponding to given 
values of the differential coefficients, the points of inflexion on a 
curve of which isolated points are given, and also in determining 
whether a series of observations is periodic. 

Example 3. — The consecutive daily observations of a function being 
0-099833, 0-208460, 0-314566, 0-416871, 0-514136, 0-605186, 0-688921, 
0*764329, show that the function is periodic and determine its period. 

From the given observations we obtain the table : — 



X 


y=/(=o) 


A 


A2 


AS 


A« 


1 


0-099833 


108627 








2 


0-208460 


106106 


-2521 


-1280 




3 


0-314566 


102305 


-3801 


-1239 


41 


4 


0-416871 


97265 


-5040 


-1175 


64 


5 


0-514136 


91050 


-6215 


-1100 


75 


6 


0-605186 


83735 


-7315 


-1012 


88 


7 


0-688921 


75408 


-8327 






8 


0-764329 
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Taking a = 3, /o =/(a) =/(3), we have, since lo = 1, 

/"(3) = «Vo-tVSVo 
= -3801-1^ -41 
= -3804 
1 .„ - 3804 



fio\ • f" (3) - 



• • /(3) ••' ^•"- 314566 
= - 0-0121. 
In the same way we find that 

(i) When tt = 4, /„=/(4) 

/"(4)=-5045, j:^/"(4)= -0-0121 
(ii) When a = 5, /„ =/(5) 

r(5)= -6-221, j^f"{h)= -0-0121 
(iii) When o = 6, /„ =/(6) 

/"(6)= -7322, ^j-/"(6)= -0-0121. 

1 

Hence in all oases 77^/"{^) = -0-0121. Putting y = f{x), this gives the 

di£ferential equation 

1 cPy 

- T^.,= -0-0121 

whence j^.^ + 0-0121 j^ = 

y = A oosO-lla' + fisinO-lla;. 
This result shows that y orf(x) is a periodic function of x, and that its 
period is 2ir/0-ll, or 57-12 days. 

Example 4. — The following table represents one of the Bessel func- 
tions, J {x) : — 



X 


y 


1-0 


0-765198 


1-1 


0-719622 


1-2 


0-671133 


1-3 


0-620086 


1-4 


0-566855 


1-5 


0-511828 


1-6 


0-455402 



dy cPy 
Find the values of -f and -7-^ for x = I -3, and hence, by substituting (with 

four-place accuracy) in the equation 

d'h/ \ dy / n^\ 

d^^+lc di + {^-l^)y = ^' 
determine the value of n. 
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Forming the table of dififerences we have 



X y 

1-0 0-765198 



A A2 ii? A* 



-45576 
1-1 0-719622 -2913 

-48489 355 

1-2 0-671133 -2558 19 

-51047 374 

1-3 0-620086 (-52139) -2184 (381) 14 

- 53231 388 

1-4 0-566855 -1796 9 

-55027 397 

1-5 0-511828 -1399 

-56426 
1-6 0-455402 

Here w = 0-1, a = 1-3, /„ =/(a) =/(l-3). 

Hence 0-1/' (1-3) = -0-052139 -J(00003S1) 
/'(1-3) = -0-52203 
•01/" (1-3) = -0-002184 -^(0000014) 
/"(I -3) = -0-2185. 
Substituting these in 

1 dy 



a'y i ay / n' \ 

d^'^ + ic d^ + v -i^)y=^' 



we have 

-0-2185 + p3(-0-5220) -f (l - jf^) (0-6201) = 

/ rfi \ 

-0-2185 - 0-4016 + yl r^i) (0-6-201) = 

-0-6201 + (l - 7^) (0-6201) = 

- 1 + 1 - (ITgji = 

n = 0. 
Hence the Bessel function in question is J^ (x). 

Example 5. — Show that the following observations are these of a periodic 
function, and find its period : — 



X (years) 


fix) 


1 


0-198669 


2 


295520 


3 


0-389418 


4 


0-479425 


5 


0-564642 


6 


0-644217 


7 


0-717356 


8 


0-783327 
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Example 6. - Observations were taken of the cooling of a boiler, and the 

results are given in the following table : - 

t (time in hours) 1-2 1'9 2-6 33 40 

B (temperature °F) 182-3 175'2 171-1 164-3 160-6 

de 
Determine as carefully as you can the value of -^ when ^ = 2-6 hours. 

MISCELLANEOUS EXAMPLES. 

1. Find a polynomial of the third degree in x which assumes the values 
37, 11, 13, 29 when x is equal to 2, 1, - 1, -2 respectively. Also find the 
minimum value of the function. 

2. A steam electric generator on three long trials, each with a different 
point of out oflf on steady load, is found to use the following amounts of steam 
per hour for the following amounts of power ; — 

Lb. of steam per hour 40-20 6650 10800 

Indicated horse-power 210 480 706 

Kilowatts produced 114 290 435 

Find the indicated horsepower and the weight of steam used per hour when 

330 kilowatts are being produced. 

3. Using first differences only, find the value of /(4817'36) if 

X 4816 4817 4818 

f{x) 693974063 694046108 694118146 
To what extent may the result be relied on ? 

4. Interpolate for a; = 35 -2967 from the following table and state the 
accuracy of the result : — 



X 


f{=c) 


35-28 


162021-40139 


35-29 


162043-22198 


35-30 


162065-04777 



5. It the sun's apparent longitude at G.M.N, on 1914 Aug. 1st was 
128° 22' 25"-8 ; on Aug. 3rd 130° 17' 15" -2 ; on Aug. 5th 13-2° 12' 7" -8 ; find 
the longitude at G.M.N, on Aug. 2nd. 

6. At the following draughts in sea water a, particular vessel has the 
following displacements : — 

Draught in feet - 15 12 9 6 

Displacement in tons 2098 1512 1018 536 

What are the probable displacements when the draughts are 11 and 13 feet 

respectively ? 

7. The following table gives the expectation of life for males of the ages 
mentioned : — 

Present age .... 5 10 15 25 35 

Expectation of life in years 50-87 47-60 43-41 35-68 2864 

Find the expectations of 3 men whose ages are 18, 20, 30 years respectively. 
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8. The H*' table of the Institute of Actuaries gives the following values 

of Ix :— 



X 


U 


10 


100000 


12 


99113 


14 


98496 


16 


97942 


18 


97245 


20 


96223 



Obtain the values of Ix when x has the values II, 13, 15, 17, 19^ 

9. The sun's apparent right ascension at mean noon is given in the follow- 



ing table :- 


- 












H. 


M. 


s. 




1914. June 1 


4 


34 


0-69 




6 


4 


54 


31-75 




11 


5 


15 


10-64 




16 


5 


35 


55-34 




21 


5 


56 


43-18 




26 


6 


17 


30-99 



Determine the apparent right ascension for 1914 June 15d. 16h. 

10. Let the number of recruits in the United States army whose height 
did not exceed x inches be denoted by y. Determine the number whose 
height exceeded 66 inches, but did not exceed 66^ inches, if 



X 




y 


62 




618 


63 




1855 


64 




3802 


65 




6821 


66 




10296 


67 




14350 


68 




17981 


69 




21114 


70 




23189 


71 




24674 


sun's apparent declination 6 


at mean noon i 


1914, May 1 


+ 14° 


54' 19"-6 


3 


+ 15° 


30' 25" -0 


5 


+ 16° 


5' 28" -4 


7 


+ 16° 


39' 27" -6 


9 


+ 17° 


12' 20" -3 


11 


+ 17° 


44' 4"-3 



Determine the hourly difference in 9 for 1914 May 5 and 7. 
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12. In the following table s denotes the space described at the time t by 
a point in a mechanism. Tabulate the values of the velocity and the acceler- 
ation, and obtain their values when t = '45. 



t (time in sees.) 


0-1 


0-2 


03 


0-4 


0-5 


0-6 


0-7 


0-8 


0-9 


1-0 


1-1 


» (space in feet) 


0-032 


0-12 


0-281 


0-525 


0-854 


1-200 


1-725 


2-229 


2-753 


3-282 


3-811 



13. There is a piece of mechanism whose weight is 200 lbs. The following 
values of 3 in feet show the distance of its centre of gravity from some point 
in its straight path at the time t seconds from some era of reckoning. Find 
its acceleration at the time t = 2°05, and the force in pounds which is giving 
this acceleration to it. 

t 2-00 2-02 2-04 2-06 2 OS 2-10 

s 0-3090 0-4931 0-6799 0-8701 1-0643 1-2631 

14. The following values of p (lbs. per sq. foot) and B°C are for saturated 

steam. Calculate with as much accuracy as the numbers will allow the value 

dp 
of ^for9=U5°0. 



e 105 


110 


115 


120 125 


p 2524 


2994 


3534 


4152 4854 


15. From the following table of the secular variations of the magnetic 


declination determine the chief period 


I:— 




Year. 






Declination. 


1540-5 






- 7° -463 


1580-5 






- 8° -500 


1620-5 






- 6° -250 


1660-5 






- 0°-500 


1700-5 






+ 8° -228 


1740-5 






+ 15° -755 


1780-5 






+ 20° -832 


1820-5 






+ 22° -380 


1860-5 






+ 19° -364 


1900-5 






+ 14° -833 



[CH. Ill 



CHAPTER III 

THE CONSTRUCTION AND USE OF 
MATHEMATICAL TABLES 

13. Interpolation by First Differences. 

All mathematical tables are alike in that they give numerical 
values of the function for certain values of the argument, which 
are so chosen that intermediate values of the function may 
be derived by interpolation. The interval of the argument 
varies in the diiferent mathematical tables, its choice in 
each case being determined by the rapidity of variation of the 
function. In most of the published elementary tables, which are 
intended for the use of students, the interval of the argument is 
so chosen that the second difference is smaller than one unit in 
the last place of decimals retained : on this account it is not 
necessary to consider second differences in the interpolation, and 
the function can be calculated for any value of the argument 
intermediate between two of the tabulated values by a simple 
proportion in the following manner : — 

Consider a function f{x) which is tabulated for the values 
Ko, x^ + w, x^ + lw, ... of the argument. Then, since the function 
increases uniformly, when digits beyond a certain place of decimals 
are neglected, we have 

/(xg + w') -/(a-'o) _ w^ 

/{Xo + W)-f{Xo) ~ w 

where (i<w'<w. 

From this it follows that 

w' 

f{x,-\rw) =f(x„) + — {f(x„ + w) -/(x,) }, 

so that, if the first difference /{x^ + w) -f{x^) is known, the value 
of the function for any value of x between »■„ and x^ + w may be 
found. This equation is known as the "Rule of Proportional 
Parts." Geometrically, what we do is to replace the curve through 
two consecutive points by the chord joining them. 
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If the first differences are only sensibly constant, the result 
obtained will be subject to a certain error. If w' = nw where 
0<n< 1, then the error may be put in the approximate form 



■ii?f"{x + n'w) 



where < w' < 1 . 



Example 1. — In the following table the entry is the cube root of the 
corresponding argument. Find the cube root uf 1619°2o. 
Argument. Entry. A 

1619 11-7421858 

24171 

1620 11-7446029 

24161 

1621 11-7470190 

If we desire accuracy up to the fifth decimal place, then we may regard the 
first difierences as constant. 

Here x-o = 1619 w = \ i«' = 0-2.5. 

Hence ^1619-25 = 11-74219 + AV (0-00242) 
= 11-74280. 

Example S. — If a; = 1000 show that, when second differences are 
neglected, the error committed in computing by this method the value of 
logjj(1000 + J«), where -« ic •« 1, will have no eflfeet on the sixth decimal 
place of the logarithm. 

Here log (a; + w) = log .t; + iw { log (a: + 1 ) - log a: } 

and the error committed is 



where 

Now 



j;=_L___^/ (a; + w,) 
•«= lo' -= 1 and /(a;) = log]oif. 



f"{x + w') 



^a logic (x + M-') 



[x+w'y- 

where m = modulus of common logarithms = 0-43429. 
But the maximum value of w{\- ro) is } since 

M(l-l«) = i -(4-!«)2<i. 

Hence the error committed is less than 

0-43429 



E 



■ (1000 + ?«')^ ' 

1 
* 16 . 10" 
«: 0-0000000625, 
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which shows that the error committed does not affect the sixth decimal place 
in the logarithm of the number (1000 + i«). 

Example S. — The differences between the six-figure logarithms, to base 10, 
of the successive integers between 5000 and 5100 are sensibly constant. 
Prove this and determine the value of the constant difference. 



14. Case of Irregular Differences. 

When the second differences cannot be neglected, however, 
this method cannot be applied. This happens frequently when 
the rate of change of the function is large compared with that of 
the argument, as is the case with the logarithmic sines and 
tangents of small angles. We must then >ise the formulae of 
Chapter II., including differences to the order required. 



Example l.—Tlo obtain log sin 0° 16' 24" '5. 
From Sohron's tables we have 



A2 A' 



ii sin 0° 16' 10" = 7 -6723450 

44543 
i sin 0° 16' 20" = 7 -6767993 -452 

44091 9 

i>sin0°16'30" = 7-6812084 -443 

4364S 
L sinO° 16' 40" =7 -6855732 

Using Gauss' formula 

J., , y .^ n(n-\) , (n+\)n(a-\) . 
/((i + n!«)=/„ + m«/i + — 2^ S^J^ + g^ asy^ 

we obtain 

£sinO° 16' 24" -5=7-6767993- 1 
19841 
56 
= 7-6787889 

This is the answer correct to 7 places, the value correct to 10 places being 
found from Vega to be 7-6787889383. 

The value obtained by the " Rule of Proportional Parts" is 7-6787834. 
This result is correct only to 4 places of decimals, and shows the necessity of 
interpolating by some such method as the above the value of the logarithmic 
sine of a small angle. This may be shown to be also true in the case of the 
logarithmic tangents of small angles. 

Example ;?. — Compute the value of log tan 0° 15' 35". 
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15. Construction of New Tables. 

In spite of the large number of different mathematical tables 
that have been constructed to assist the computer, further 
research in the various branches of mathematics and the kindred 
sciences shows that these are not sufficient for all purposes. A 
research student often requires new tables for the question which 
he is discussing, and these he must construct for himself. 

The first question to be settled is the degree of accuracy to 
which the tables are to be computed, i.e., the number of decimal 
places to be retained. This varies in different cases ; but it is 
always desirable to extend the original calculations beyond the 
number of places which are ultimately to be retained, in order to 
avoid the cumulative effect of the neglected digits. 

The general plan of the calculation is to compute first with great 
accuracy the values of the function for a series of widely -spaced 
values of the argument. For this purpose infinite series are generally 
employed. The values of the function for a more closely-spaced 
series of values of the argument are then derived from these by 
means of the formulae of Chapter II. : this second process is called 
Suhtabulation. We shall consider these processes in §§ 17, 18. 

16. Tabulation of a Polynomial. 

The simplest ease to be considered is the tabulation of a polynomial. 

If the degree of f(x) is n, it will not be necessary to compute its value 
directly for more than {n + \) values of the variable, as may be seen from the 
following examples. 

ExanvpU l.—hei f{x) = ofi + Zx^ + 2x-\. 

When a;=0, I, 2, 3, we have f(x)= - 1, 5, 23, 59. We thus obtain the 
following table :— 

X f(x) A A2 A3 



1 5 12 

18 6 

2 23 18 

36 

3 59 

Since f{x) is of the third degree in x, its third differences will be constant ; 
in this case they will be equal to 6. Extending the column A' by inserting 
6*8 we can then extend the column A^ by simple additions or subtractions. 
A repetition of the same process will enable us to extend the column A, and 
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thereafter to dcterinine the values of f(x) for all integral values of x. The 
resulting table is 



X 


/(^•) 


A 


^^ 


A3 


-3 


-7 


6 


-12 


6 


-2 


-1 





-6 


6 


-1 


-1 








6 





-1 


6 


6 


6 


1 


5 


18 


12 


6 


2 


23 


36 


18 


6 


3 


59 


60 


24 


6 


4 


119 


90 


30 


6 


5 


209 




36 





By this simple method we are enabled to determine the values of a poly- 
nomial f{x) for all positive and negative integral values of x. We can then 
use any of the formulae of interpolation already established to obtain the 
values of /(x) for fractional values of x. The method may be extended to 
functions which are not rational or integral, but in this case it will be 
necessary to make the interval between successive values of x so small that 
the diflferenoes of some definite order are constant throughout the range of 
values under consideration. 

Example 2. — Form the table of values of /(.x) = x^ 4- 03? - Ix"- + x - 4 for the 
values 0, 1, 2, 3, 4 of the argument .f, and then deduce the values of /(x) 
corresponding to the values - 1 , - 2, - 3, - 4 of x. 

17. Calculation of the Fundamental Values for a 
Table of Logarithms. 

We shall now illustrate the process of computing tables of 
functions other than polynomials by describing the calculation 
of a table of logarithms. 

The first stage in this consists in determining the logarithms 
of all the prime numbers between certain limits. It is evident 
that when the logarithms of the primes are known, the logarithms 
of all composite numbers can be derived from them by mere use of 
the formula 

log a 6 = 102 a + log h. 
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In order to obtain the logarithms of the primes we take two 
equations whose roots are integers and which differ only in their 
constant term. Such a pair are 

ar-3a3 + 2 = 
a!=-3a;-2 = 0. 
Now a;'-3a! + 2 = (a;-l)=(a; + 2) 

Hence 

log (as - 1) + i log (a'+ 2) - log (a;+ 1) - ^ log (a; - 2) 

«■' - 3 a: + 2 



= iloo 



1 + 



■3a;-2 



, , X' — ^ X 

= i log 

1 



a? -ix 



which, by use of the ordinary logarithmic series, may be written, 
if the logarithms are all to the base e. 






2 
a? -Zx 
This is called Borda's Formula. 

Putting a; = 5, 6, 7, 8, we obtain the four equations 

log 2 - f log 3 +i log 7 = 0-0181838220854376... 

I log 2 + logs- log 7 = 0-0101013536587597... 

- 2 log 2 + 2 log 3 - i log 5 =0 0062 11 2599992784... 

- f log 3 + ^ log 5 + log 7 =0-0040983836020895... 

Solving these four equations we obtain the logarithms (to the 
base e) of the prime numbers 2, 3, 5, 7, viz. ; — 

log 2 = 0-6931471805599453... 
log 3 = 1-0986122886681096... 
log5= 1-6094379124341003... 
log 7 = 1-9459101490553133... 

Further substitution for x in Borda's formula will give immedi- 
ately the logarithms of the other prime numbers. For example, on 
putting a; = 9 the left-hand side becomes 

Iog8 + |logll-logl0-^log7. 
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Of these log 7 has been found ; 

logs =3 log 2 = 2-0794415416798359... ; 
and log 10 = log 2 + log 5 = 2-3025850929940456.... 

Hence we obtain log II = 2-3978952727983705 .... 

It is to be noted that these logarithms are to the base e. But 
as we have now found log^ 10, we can at once derive its reciprocal 
logioe = 0-4342944819032518..., and all the above logarithms 
can now be converted to the base 10 by multiplying them by this 
" modulus." 

We thus obtain the table 

logio 2 = 0-3010299956639812 
logio 3=0-4771212547196624 
logio 5=0-0989700043360188 
log,o 7 = 0-8450980400142568 
logioll = 1-0413926851582250 
etc. 

For the calculation of the logarithms of the larger primes Haros' Formula 
will be found more useful. In this case consider the equations 

a;*-25a;=' + 144=0. 
From these we deduce 

'°S(a; + 3)(x-3) (a; + 4)(rE-4) ^^°^ x^-2bx^+lii ^ ^°^ {x*- 25x^+72) +72 
whence, if the logarithms are to the base 10, we have 

Iog(a: + 5) = {log(a; + 3) + log(j,--3) + log(.-c + 4) + log(a;-4)}--[Iog(a;-5) + 21oga:} 
r 72 / 72 \3 \ 

"^'"\.-«r'-25a;2 + 72 + *W-23a;= + 7i/ "*"■■ I 
where ni=logiQe. 

Example 1. — To obtain the logarithm of the prime number 43 to the 
base 10. 

Let a; + 5=43. 

Then .-k = 38 

log 43 = log 4 1 + log 35 + log 42 + log 31 - log 33 - 2 log 38 
- 0-8685889638065036 {0-0000351372402040 + -0000000000000145}. 
But the logarithms of 41, 35, 42, 34, 33, 38 are assumed to have been 
previously calculated. Hence 

log 43 = 1 -6127838567197355 - 1 -5185139398778875 
1 -5440680443502756 3 1595671932336203 
1 -6232492903979005 0-0000305198190724 
1-5314789170422551 
= 1-6334684555795865 
Example 2.— Show that logn, 131 = 2-1172712956557643. 
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18. Amplification of the Table by Subtabulation. 

Having now obtained the logarithms of the prime numbers 
within a certain range, we derive from them, by mere addition, the 
logarithms of the composite numbers which lie in this range, and 
thus obtain the logarithms of all the integers within the range. 
We shall now show how, by subtabulation, the extent of this table 
may be increased manyfold. It will be supposed that 7-place 
accuracy is required. 

For example, suppose the logarithms of the integers 31, 32, 33, 
34, 35, 36 have been obtained directly ; we shall show how the 
logarithms of all the integers between 330 and 340 can be derived 
from them. 

Neglecting the characteristics of the logarithms, which can 
always be written down by inspection, we form the following 
table : — 



X 


log X 


A 


A'-' 


A3 


A* 


31 


■49136169 


1378829 








32 


■50514998 


1336396 


-42433 


2535 




33 


■51851394 




- 39898 




-223 






1296498 




2312 


( - 208) 


34 


•53147892 


1258912 


- 37586 


2120 


-192 


35 


•54406804 


1223446 


- 35466 






36 


•55630250 











Putting a = 33, n = -^^, io=\ in Gauss' formula 

71 ( 1 -«.)„„ . W ( 1 - «■•') 

/{a + nw) =/o + n S/, ;^-j — 6-/o " 



3! 



871 



^^ISLZ^^^.BV, 



we obtain 

/(33-1) =/; + 0-1 Vj - 0-045 6% - 0-0165 37^ + 0-0081 fi h*J\ 
= •51851394-38 
129650 1 
1795 
= ■51982800. 
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Thus to 7 decimal places log 331 = 2-5198280. But it is not 
necessary to calculate the values for all the numbers between 
330 and 340 in this way. It is, in fact, simpler to calculate 
log 332 from the formula 

f{a + -^ w) =/{a + -^V «>) + {/(« + tV "') -/(«) } + S= 
where S., =f{a +-^w)- If {a + ^w) +/{a) 

= 0-01 sv; + 0001 s'/j - 0-001 /^ sv' 

and then to calculate log 333 from the formula 

/(« + TU W') =/(« + To "') + {/(« + TO »'') -/(« + tV «>)} -f Ss 
where 83 =/(a + f ^ w) - 2/(a + ^-^ w) +/(a + ^^ w) 

= 0-01 sy„ + 0-002 8v^ - 00017 jj. sy^ . 

The complete table of these quantities S is as follows : — 

s„ =0-01 8y„ + 0-001 sy^ - o-ooio /* sv; 
S3 =0-01 gy(,+o-oo2 sy^- 0-0017 ^^sy 
8, = 0-01 sy„ + 0-003 sy; - o-ooi9 /^ sy 
8, =0-01 sy + 0-004 8y - 0-0021 /. 8yJ 
8, = 0-01 8y + 0-005 sy - 0-0020 /* sy 
8, = 0-01 sy + 0-006 sy - 0021 ^a sy; 
Sg = 0-01 sy + 0-007 sy^ - 0-0018 /* sy; 
8s = 0-01 sy + 0-008 sy^ - 0-0017 yx sy^ 
8i„ =0-01 sy + 0-009 sv; - o-ooi 2 /x sy^ 

The numerical work is arranged as in the following table, in which 
the logarithms required are found in the last column. 



— i 


J9898-0 + 2 


312-0-2 


08-0 




129649-8 


2-51851394 


tV + 


1795-4- 


38-1- 


1-7 = 


-1755-6 


131405-4 


2-519827994 


1% ~ 


399-0 + 


2-3 + 


0-2 = 


- 396-5 


131008-9 


2-521138083 


3 

To- ~ 


399-0 + 


4-6 + 


0-4 = 


- 394-0 


130614-9 


2-522444232 


4 
To - 


399-0 + 


6-9 + 


0-4 = 


- 391-7 


130223-2 


2-523746464 


6 

To ~ 


399-0 + 


9-2 + 


0-4 = 


- 389-4 


129833-8 


2-525044802 


c 
TT) - 


399-0 + 


11-6 + 


0-4 = 


- 387-0 


129446-8 


2-526339270 


7 
TT7 ~ 


399-0 + 


13-9 + 


0-4 = 


- 384-7 


129062-1 


2-527629891 


8 

TU - 


399-0 + 


16-2 + 


0-4 = 


- 382-4 


128679-7 


2-528916688 


9 

T(r ~ 


399-0 + 


18-5 + 


0-4 = 


- 380-1 


128299-6 


2-530199684 


10 


399-0 + 


20-8 + 


0-2 = 


- 378-0 


127921-6 


2-531478900 
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The method followed in eonstruoting this table will best be shown by 
considering two examples. 

In the line ^ 

+ 1795-4 is the numerical value obtained for 0-045 5'/o 

- 381 „ ,, „ -0-0165 «s/j 

- 1'7 ,. „ „ 0-008 M«V^ 
1755-6 is the sum of the three preceding numbers. 

129649-8 is the numerical value obtained for 001 Sf,. 
131405-4 is the sum of the two preceding numbers, and therefore 
represents 

0015/^+0-045«V„-00165 5Vj + 0008m«Vj- 
2-51851394 is the numerical value of /„ 

2-519S27994, which is obtained as the sum of the two preceding 
numbers, represents 

/, + 0-01 S/^ + 0-045 5V„-0-0165S3yj+0-008/x5V, . 
In the line ^ 

- 399-0 is the numerical value obtained for O'Ol S^/o 
+ 2-3 „ „ „ 0-001 SV^ 

+ 0-2 „ „ „ -0-001 /iS-i/j 

-396-5, which is obtained by summing the three preceding numbers, 
represents 

0-01 SPfo + 0-001 «V^ - 0-001 M S-*/^ or Sj. 

131405-4, which has previously been obtained, represents 

131008-9, which is obtained by summing the numbers immediately 
above it and to the left of it, represents 
/(a + Aw)-/(a) + «2. 
Now 2-519827794, as we have already seen, represents f{a + -^w). 

2-521138083, which is obtained by adding the numbers immediately 
above it and to the left of it, represents 

The results in the last column are the logarithms of the numbers 
330, 331, ... 340. They cannot, of course, be relied on up to the 
last figure, but -will be found correct up to seven decimal places. 
This accuracy may be obtained by taking one decimal more in the 
first difference than in the tabular function, one decimal more in 
the second difference than in the first, and so on, when these can 
be obtained. A check on the accuracy is obtained by carrying 
on the process until log 340 is computed. Comparison of this 
result with the value for log 34 will show whether the computa- 
tion has been accurately performed. 
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We may remark that it is an improvement on the foregoing 
process if the diiferences required for the subtabulation are 
computed independently from the expressions in series, instead of 
being derived by forming a table of differences. 

Example 7.— Given log 24 = 1 -3802112 

log 25 = 1-3979400 
log 26 = 1-4149733 
log 27 = 1 4313638 
log 28 = 1-4471580 
log 29 = 1-4623980 
obtain the logarithms of the numbers 261, 262, ... 269. 

19. Radix Method of Calculating Logarithms- 

When the logarithm of a number is required to a large number 
of decimal places, its determination either by means of an infinite 
series or by interpolation becomes very laborious and liable to 
error. A much simpler and more useful method in such a case is 
that known as the Radix Method. 

Let N be the given number whose logarithm is required to n 
places of decimals. Then since 

where N^ is a decimal, x an integer, and a = 0, 1, 2, ... 9, 

logio iV = a; + log,o a + logi„ {\+N„). 

Thus, to determine logmiT, it will be sufficient to compute logjoO 
and logio (1 + N^. The former has been obtained above. We shall 
now show how to evaluate the latter. 

Let r^ denote the decimal consisting of the first m figures in 
Nf„ and let N-^ be the difference of N^ and r,. Then 

(l+J7„)-(l+rO = iV,. 

Now divide N-^ by (l+»-i) until there are m figures in the 
quotient r^. If we call the remainder iVj, then 

N,-{\+r,)r„ = N^. 
Similarly, if the remainder, after dividing N^ by (1 H-**]) (1 H-rj) 
until there are m figures in the quotient 7-3, is denoted by N,, we 
have 

and in general, 
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These equations may be written in the form 
(l+If,) -iV, =(l+n) 

iY, -W, =(l+r.)(l+r,)-(l + rO 

-fl^,-2-^p-, = (l +n) (1 +r-,) ... (1 +r,_0 

-(l+r,)(]+r,)...(l+r,_,) 

-^P-i--^,. =(l+r,)(l+r,)...(l+rj 

-(l+rO(l+r,)...(l+r,_,). 
Adding, we obtain 

(l+i7„) -N^ =(l+r,)(l+r,)...(l+r,) 

whence 

1 +iro =(1 +n) (1 +r,) ... (1 +r,) + #^. 

If, therefore, the remainder JFj, has no significant figures in its 
first n decimal places, 

log(l+i\r^) =log(l+r,) + log(l+r2)+ ... +]og(l+r,,). 

The quantities (1+r,), (l+r^), ... are called the radices. As 

they will be of the form ( 1 + Yrv) ''^liere k and I are integers, and 

as k is always small compared with 10' we may compute log (1 +rp) 
by means of the logarithmic expansion 

log(l+»-,) = r,-^r/ + ^r/-iV+... 

The values of the logarithms of these radices when A; = 0, 1, 2, ... 
999, and 1 = 3, 6, 9, 12, 15, 18, 21, 24, as well as of the numbers 
1, 2, ... 9, have been computed to twenty-four places of decimals.* 

The inverse process of determining the antilogarithm corre- 
sponding to a given logarithm is equally simple. When the 
logarithm is not in the table, we have merely to take out the 
next lower and subtract it from the given one. Then take out 
the next lower and subtract it from this remainder; and so on 
until the remainder consists of zeros. The product of the radices 
corresponding to the logarithms taken out is the number required. 

* Gray, Tables for the formation of Logarithms and Antilogarithms to 
twenty-four or any less number of places. London, C. & B. Layton (1900). 
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We shall illustrate these two methods in the following 
examples : — 

Example 1. — To determine logio43 to 24 places of decimals. 

Here a; = l, a =3. 

Taking to = 3, as in Gray's tables, we have 

1 + i^o = 1 • 433 333 333 333 333 333 333 333 
l+ri=l-433 
iVi = ■ 000 :<33 333 333 333 333 333 333 
r2=0- 000232 

i^2=0- 000 000 877 333 333 333 333 333 
(l + ri)r-2= 0-000 332 456 i.e., 0-000333 333-0' 000 000 877 
(l + ri)(l+r2i = (l + !-,) + {! +?-i)r-, 
= 1-433332456. 
The above will be sufficient to indicate the general theory. The process 
can, however, be simplified by arranging the work in the following manner : — 
1 - 433 333 333 333 333 333 333 333 
1-433 

Di£F. =0-000 333 333 333 333 333 333 333 

1 - 433 ) 0- 000 333 333 333 333 333 333 333 ( 0- 000 232 
286 6 



46 73 
42 99 


3 743 

2 866 



1-433333333 

0-000000877 



877 = Remainder. 



Uifif. = 1-433 332 456 = New Divisor. 

1 - 433 332 456 ) - 000 000 877 333 333 333 333 333 ( 0- (000)= 612 * 
859 999 473 6 



17 333 859 73 
14 333 324 56 

3000535173 
2 866664912 



133870261 = Remainder. 
1-433333333333333333 

0-000000000133 870261 
DifiF. = 1-433 333 333 199 463 072 = Ne w Divisor. 

Should the original dividend not be known beyond the 24th place of 
decimals, we should now have to proceed by a method of contracted division. 
In the present case we know the figures after the 24th decimal, and may 
therefore continue as above. But we shall adopt the contracted method. 

♦ 0-(000)2612 is written in place of 0000000612 
0- (000)' 093 „ ,, 0-000000000093 and so on. 
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The new divisor oontaina 19 Bgures ; the new dividend 
•tOOO)' 133 870261 333 333 
contains only 15 significant figures after the point. We must therefore cut 
off 4 figures from the divisor. Thus we obtain 

1 • 433 333 333 199,4,6 ) 0' (000)' 133 870 261 333 333 ( 0- (000f093 * 

128 999999987951 

4870261345 382 
4299999999598 



570261 345784 = Remainder. 
1-433 333 333 333 333 333 333 333 
0- 000 000 000 000 570 261 345 784 



Diff. = 1-433333333332763071 987 649 = New Divisor. 
As the new dividend 0- (000)^570261 345784 contains only 12 significant 
figures and the divisor 25 we must out ofl' 13 figures from the latter. It 
then becomes 1-43333333333. As these are the figures in the original 
dividend there will be no further new divisors. This will usually happen 
when the fourth division is reached. We may therefore determine the 
remaining figures by continued contracted division. 

1 - 433 333 3S3 33 ) 0- (000)^ 570 261 345 784 ( - (000)^ 397 856 752 873 
430000000000 



140261345 784 
129 000000000 

11261345 784 
10033 333333 

12-28 012451 
1146666666 

81345 785 
71666 667 



9 679118 
8 600000 

1079118 
1003333 

75785 
71667 

Tii8 

2867 

1251 
1146 

105 
100 

5 
4 
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The remaining radices are thus 

!• (000)'' 397 
1- (000)5 856 
1- (000)8 752 
1- (000)' 873 

Using Gray's tables of the radices we have 

logs =0-477121254719662437295028 

log 1-433 =0- 156 246 190 397 344 475 994 693 

log 1-000232 = 100744 633875 845680796 

log 1 - (000)2 612 = 265 788 141 593 627 087 

log 1 - (000)3 093 = 40 389 386 815 124 

log 1 • (000)^ 397 = 172 414 909 316 

log 1 - (000)5 856 = 371 756 077 

log 1- (000)6 752= 326589 

log 1- (000)' 873= 379 

Sum = 0-633468 455579 586526405 089 



Hence since x = l the required logarithm is 

log 43 = 1 - 633 468 455 579 586 526 405 089. 

Example 3. — To determine the antilogarithm of 

1 • 120 619 882 724 133 396 646 450. 
As the integral part only a£Fects the characteristic we shall neglect it for 
the present. 

Using Gray's tables we then obtain 






120619 882 724133396 646450 


log 1-320 =0-120573931205 849 868 472706 


Difif. 


45 951518 283528173 744 


log 1-000105 = 


45 598526719080137349 


Diff. 


352 991564448 036395 


log 1 -(000)2 812 = 


352 646976130787551 


DifiF. 


344588 317248844 


log 1 -(000)3 793 = 


344 395 524 012726 


Dift 


192793236118 


log 1 -(000)^43 = 


192 392455483 


Diff. 


400780635 


log 1 -(000)5 922 = 


400419512 


Diff. 


361 123 


log 1 -(000)6 831 ^ 


360899 


Diff. 


224 


log 1(000)' 516 = 


224 


Diff. 






Multiplying these radices together we have 
1 -(000)^ 443 X 1 -(000)5 922 x 1 -(000)6 331 x 1 -(000)' 516 = 1 -(000)^ 443 922 831 516. 
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This product of the last four factors can always be written down by 
inspection. The other multiplications may be arranged as follows :— 
1 ■ 000 000 000 000 443 922 831 516 

1 • 000 000 000 793 

1 • 000 000 000 000 443 922 831 516 

700000000000310 

90000000 000040 

3000000000000 

1 • 000 000 000 793 443 922 831 866 
1-000000812 



1 • 000 000 000 793 443 922 831 866 

800000000634 755138 

10 COO 000 007 934 439 

2000000001586 888 

1 ■ 000 000 812 793 444 567 108 331 
1- 000 105 

1-000000 812 793444567 108 331 

100000081 279 344 456 711 

5 000004 063967 222836 

1 - 000 105 812 878 787 878 787 878 
1-320 

1 • 000 105 812 878 787 878 787 878 

300 031 743 863 636 363 636 363 

20002116 257 575757575 757 

1 ■ 320 139 672 999 999 999 999 998 



Since the characteristic in the given logarithm is 1, there will be 2 figures 
before the point in the corresponding antilogarithm. Thus the required 
result is 

13- 201 396 729 999 999 999 999 98. 

Example 5.— Given 7r = 3- 141 592653589793238462643 show that 
log 7r = 0-497 149 872694 133854351268. 

Example i. — Show that the antilogarithm of 
0-434294481903251827 651129 is 2-718 2818-28459045235 360290. 

20- Inverse Interpolation. 

Suppose the values of a function have been tabulated, and that 
it is required to find the value of the argument which corresponds 
to some given value of the function, intermediate between two of 
the tabulated values. The process by which we obtain this value, 
which may be called the "antif unction," is known as Inverse 
Interpolation. 
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liet/ia + nw) be the given value of the function. Then 

f{a + n w) =/„ + n 8/^ + -^-^y ^ Sy, + . . . 
n{n-\)...{n-r+l) 

' • 2 

In this equation f(a + m w) is given, and as the values of 
fojfi, ... have been tabulated, the values of S/^, 8^,, ... are also 
known. Hence the only unknown quantity in the equation is n, 
and this is the number required for the determination of the anti- 
function (a + WM)). 

To obtain it we must solve this equation. It has been pointed 
out already that the differences of the function under consideration 
usually become sensibly constant at some definite order, say r. 
The higher differences will therefore be zero, and this equation then 
becomes one of degree r in n, It can, therefore, be solved by a 
method of successive approximation, as follows : — 

If we write the above equation in the form 

n = l=fi 

and neglect differences of the second and higher orders, we get as 
a first approximation the value Wj where 

Let us next take second differences into account and put n-^ for n 
in the denominator of the expression on the right-hand side. Then 
the next approximation n.^ is given by 

/-/o 

When third differences are taken into account, n^ is substituted in 
the denominator for TC], giving 

n = — f-fo 

' 84 + H'*^ - I) Sy, + i{n,-l) (n, - 2) Sy " 

Proceeding in this way we obtain closer and closer approxima- 
tions to the value of n. The method, though rather laborious, has 
the advantage that any error introduced in the computation of any 
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one approximation will not cause the final result to be wrong, 
being itself rectified in the higher approximations. Thus the sole 
effect of the eri-or introduced is merely to make the approximation 
less rapid, and therefore to increase the amount of labour involved 
in the accurate determination of n. On this account the process is 
a safe one to employ. 

Other formulae for determining the antifunction may be 
obtained by reverting any of the formulae of interpolation. We 
shall exemplify this by reverting Stirling's formula 

Neglecting differences of higher order than the fifth, and arranging 
in powers of n, we have 

/(a + nw) =/„ + {f.S/„-^f. SV; + ^fi 8y„) n + (J Sy„ - ^r «'/») n' 
+ {ll^ 8y„ - J^ /. Sy ) n' + Jf 8Vo n' + jh /* ^Vo n' 
which may be written 

F^An + Bn'^ G'n? + Dn' + En^ 
where 

/'=/(a + nM))-/o 
^=/-'S/„-i/x8y„ + .V/^Sy„ 

^=isy„-.v8y„ 
c = |/.^8y„-J^/.8y, 

Reverting this series we obtain successively 

F 

(i) -=x 

F B (F\- 

F B (F B fFVV C /Fy 

(iii) «= y - -J \T - T U; / " :i U z* 

F B /FWi „/BV C\(Fy 
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<">•.=^f[^4(^^{<4)'^}{zn' 

~ :a XT ~ 1\a) I a\a) 

-{»¥-^-=(t)'}©' 

W"-7-7(j)nH7)-7Mxl 

A more convenient form * may be obtained by putting 

— = 71, and — = r and rearranging the terms. 
A A 

We then obtain the formula 

n = wj +/i Wi +/, »i/ 4-/3 Hi^ +/4 n/ 

where 

/i = - 5r + 2 (5r)= - 5 (5r)= + 1 4 (5r)^ 

/2 = Cr{-l + 5^r-21(5r)=} 

J, = Dr{-l+6Br) + 3{Crf 

f,= -Er. 

If the first derivate of the function f{a + nw) is tabulated, the 
values of the quantities A, B, C, D, E may be obtained by forming 
the successive differences of wf {a + n w). As fourth dififerences of 
the derivate are of the same order as fifth difi'erenees of the 
function we need not extend the table beyond this order of 
differences. 

* Van Orstrand, Phil. Mag. (6) 15 (1908), p. 630. 
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Argument. Entry. A A^ A' A* 

a — 3w wf'_.i 

a-w ■*«/'_] wS-y'_] w^*f-i 

a + 3 w wy'j 
Differentiating the equation 

/{a + nw) =/„ + ilw + Bn' + Gv? + 2>n^ + Ew' 
we see that 

A=wf\; B=^\w'fr, G=^\wT; D = ^-,v^fi'^ 

Now from Stirling's formula we have 



2 / 2 "1 \ 

w/'(a + ww) = w/'o + ni«jii8/'„+ '%^wK-f\^ "31" '^V-^^'f'" 

+ 4, "'g/o 
w=/" (a + w w) = i«/iS/'„ + n?«8y'o +■ ^g~ w /iSy'o + "^^^ !^8y'o 

i(fif^"{a + nw)= wSy'o 

Putting w = in these equations we obtain 
A=tv/'{a) 

B = J wT (a)^iwi^ s/'o - A «> /* sy '„ 
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A third method of determining the antifunction is to determine 
rii as in the first method, and then, using shorter intervals, to 
retabulate the function for values of oi in the neighbourhood of rii. 
The second approximation is now obtained from this new table in 
the same way as n^ was obtained from the original table. One or 
two repetitions of this process will give the value of w to a very 
high degree of accuracy. 

An approximation to the value of the error committed by computing n 
by the first method will next be obtained. 

Since 



sf, +- 



M-1 

f-fo 



5/ +^ii^/"(^ + i'w) 



■^{- 



n-1 

2! 



^f"{x + n'w) 



«/i 



r 



the numerical value of the error committed is 

/-/o 



^1 = 



«4 ■ 

ii{n-l) 
2l 



n-1 f"(x + n'w) 

"Tr'" ■ 



II? 



*4 

f"(x + n'w) 
Ji~fo 



Example 1. — The following table gives the logarithms of the distance 
of Venus from the earth. To find when the logarithm had the value 
9-9351799. 



Date. 




log. 


A 


A2 


A3 


A^ 


A» 


1914 Aug. 


18 
20 


9 '9617724 
9-9543645 


-74079 
-75576 


-1497 


-46 








22 


9-9468069 


-77119 


-1543 


-43 


3 


1 




24 


9-9390950 


(-77912) 
- 78705 


-1586 


(-41) 
-39 


4 


(1) 





26 


9-9312245 


-80330 


-1625 


-35 


4 






28 


9-9231915 


- 81990 


-1660 









30 9-9149925 
We shall use the second method described above. 
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Let 


/o= 9-9390950. 






Then 


^= -0-0077912 + 0-0000007= -0-0077905 






£=-0-0000793 








0= -0-0000007 








D= 








A'= 








9-9351799-9-9390950 






"i- -0-0077905 


= 0-502535 






r= -64-5062 








Br= 0-0051153 








Cr= 0-0000452 








Dr= 








Er= 








./i=- 0-0050637 


/i?H =-0-002545 






'/2=- 0-0000440 


/2<= -0-000011 






/3= 


«i= 0-502535 






/.= 


,<, = 0-499979. 




Hence the time 


required is 








1914 Aug. (24 d. + 0-499979 x 48 h.) 






or 1914 Aug. 24 d. 23 h. 


59 min. 56 s. 




Example S. 


—To show that when a:= 


1000 and !« -=: 1 the error in com 


puting w from the formula 







log {x + w) - log X 
~log (a;+ 1) -logo; 
is less than 0-000125. 

Here 

I wll -w) f" {x + w') 



2 ! 
k{\- 



■/o 



Id) 



where • 
1 



2 (.T + i'/f log (a: + 1)- log X 

But by Taylor's theorem 

m 
log(a; + l)-loga:= ^. + ,;," whereO-sif" 

t«(l-«») x-\-vi" 



\E\ 



(x-^-w'f 



<\^ 



TTrtnr 



< 0-000125. 


g table 


gives the va 


X 


Ux) 


1-1 


0-47090 


1-2 


0-49829 


1-3 


52202 


1-4 


0-54195 


1-5 


0-55794 


1-6 


0-56990 


1-7 


0-57777 



Find the value of x corresponding to the value 0-55302 of J-^{x). 
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21. Various Applications of Inverse Interpolation. 

The method of inverse interpolation may also be employed to 
obtain the maximum and minimum values of a function which is 
given as a series of observations, or as a graph. In this case we 
should differentiate one of the formulae of interpolation, say, 

/{a + nw)^J,^n 8/^+ ^^^^^^^ 8y, + ... 

Differentiation with respect to n gives 

For maxima or minima we have 

/' {a + nw) = 0, 
so that n is given by the equation 

which is solved in the same way as the corresponding equation in 
last section. 

When it is desired to obtain the value of n, for which the 
gradient of the function has a particular value c, say, the equation 
for determining n becomes 

The same method may be employed to obtain the solution of 
any equation involving one unknown, whether it be algebraic or 
transcendental. For we may tabulate the values of the expression 
occurring in the equation and obtain their successive differences. 
The solution of the equation is then the same as finding the values 
of the argument corresponding to zero values of the function. 

Example i.— Approximate to the roots of the equation 
a;3 + 4a;2 + 3a;-2=0. 
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A A2 A' 

1-875 

-1-750 
0-125 0-75 

-1-000 

■ 0-875 0-75 

-0-250 

■ 1-125 0-75 

0-500 
- 0-625 0-75 

1-250 
0-625 0-75 

2-000 
2-625 0-75 

2-750 
5-375 0-75 

3-500 
8-875 0-75 

4-250 
13-125 0-75 

5-000 
18-125 0-75 

5-750 
23-875 



The positions of the roots will be found by examining the column f{x). 
We then see that when /(a;) = 0, the value of x lies between and 0-5. It 
will also be seen that f{x) is zero when x= -2, and that there is a 
possibility of another root in the interval -2-5< a- < - 1-5. We shall 
consider the various approximations to the root for which 0< a; < 0-5. 

Using the formulae of the first method of § 20 we have 

/-/o 

0-(-2) 
~ 2-625 
= 0-7619 

f-A 



'irst form the table. 


X 


f{x) 


-3 


- 2 


-2-5 


- 0-1-25 


_2 





-1-5 


- 0-875 


-1 


- 2 


-0-5 


- 2-625 





- 2 


0-5 


0-625 


1-0 


6 


1-5 


14-875 


2-0 


28 


2-5 


46-125 


3-0 


70 



2 

2-625-0-119x2-75 
: 0-8704 
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/— / 
"■^ = «/,+4(n2-l)S'/i + iK-l)(n2-2)«»/, 

2 

= 2-625-0-0648x2-75 + 0-0244x0-75 

= 0-8113 

/-/o 



2 

= 2-625- 0-0944 X 2-75 + 0-0374 X 0-75 
= 0-8356. 

The first four approximations to the value of n are therefore 
n, = 0-7619 
n., = 0-8704 
Us = 0-8113 
n, = 0-8356 

But the tabular interval is one-half. Hence the first four approximations to 
this root are 

Xi = 0-3810 

^2 = 0-4352 

X3 = 0-4057 

x-4 = 0-4178 

•The correct value of the root is 0-4142. Subtracting each of the values 
obtained from 0-4142 vi^e get 

0-4142 -X, = +0-0332 
0-4142 -.r^ = -0-0210 
0-4142 -a;3= +0-0085 
0-4142 -:r4 = -00036 

These differences show that each step gives a much closer approximation 
to the root, and that by continuing the process we can get a result which 
shall be correct to any given number of decimal places. 

Examples of the second method, in which only first differences are taken 
into account, will be found in any text-book which discusses the solution of 
numerical equations. 

Example S. — Approximate to the real roots of the equations 
(1) e-'-x^=0 

(ii) e^ + 3i?-4 = Q 
(iii) 0-5a-i-=- 121ogi(,a; + 2sin2a,-=0-921. 
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MISCELLANEOUS EXAMPLES. 

1. If y=a;tin-'a;, calculate a table of first differences near a. = 0'43 for 
increments 0'0002 in x. 

2. The following table gives the moon's right ascension 



Date. 






R. 


A. 


1914, March 5d. 


Oh. 


5h. 


3 m. 


7 -26 8. 




10 li. 




26 m. 


58 '968. 




20 h. 




51m. 


14-83 s. 


6d. 


6h. 


6h. 


15m. 


48-968. 




16 h. 




40 m. 


34-76 s. 


7d. 


2h. 


7h. 


5 m. 


25-36 s. 



Complete the table for every hour from March 5d. 20 h. to March 6d. 6h. 

3. Prove that in a table of logarithmic tangents to base 10, the difference 
for one minute in the neighbourhood of 60° is 0-00029. 

4. Determine the values of (1) log tan 1° 29' 33" 

(2) log sin 1° 15' 12" 

5. Show that logn, 1031 = 3- 013 258 665 283 516 546 909 664. 

6. Show that the number whose logarithm to 20 places of decimals is 
given by 0- 004 892 890 303 239 039 78 is 1 -01 133. 

7. Prom the following table determine the value of oc for which f{x) has 
the value 0-22389. 

X f(x) 

1-5 0-51183 

1-7 0-39798 

1-9 0-28182 

2-1 0-16661 

2-3 0-05554 

8. Approximate to the real roots of the equations 



(i) 


xe"- + 2a;-5=0. 




(ii) 


ex-e-» + 0-4a;-10=0 




(iii) 


2a:''i-3a;-16 = 




(iv) 


2-42 a^- 3-15 log^x- 20-5=0. 




9. From the series 






hovx = 1 - 22_43 


^8 3.12 




"*" 22.42.6=. 8^ ~ 2^.42.6^ 8=. 


10M2= 



compute a table of ber x from a; = 0-0 to a; =5-0 at intervals of 0-5, correctly 
to 6 places of decimals. 

Hence, by subtabulation, compute a table of ber a; from a:=2-0 to z=3-0 
at intervals of O'l. 



62 INTERPOLATION AND NUMERICAL INTEGRATION [ch. hi 



10. The confluent hypergeometric function W^ ^^ {z) is computed for large 
values of z from its asymptotic expansion 

^ {rrfi -{k- if) K - (fc - ^)'} .■.{m^-{k-n + m \ 
n ! a" I 






while for small values of z it is computed from the formula 

T(-2m,) T{2m) „ 

"'.. .» (^) = T{i-m-k) ^*. » (^) + r(i + m-ft) ^.. -... (-^' 

where 

C i + } n-k {i + m-k){^ + m-k) ^ 

M^, „. (z) = 2^+"- e-i'^ll + 1! (2m + l) ^ + 2! (2m + l)(2m + 2) ^ + •••/ _ 

Compute tables and draw graphs of W^. ^Jfi) for positive values of z in the 
cases 

i=-3-l, m=+2-2 

i=-0-l, m=+0-2 

i=+2-2, m=+0-2. 



CHAPTER IV 

NUMERICAL INTEGRATION 

22. Introduction. 

The method of evaluating the definite integral of a function 
from a series of numerical values of the function is called 
Mechanical Quadrature or Numerical Integration. 

If a function is given by its graph, the simplest way of deter- 
mining an area bounded by the curve, two given ordinates and a 
given abscissa, is by means of a planiraeter. This method would 
naturally be used in the case of the indicator diagrams of steam, 
gas, or oil engines, or the cards from certain hydraulic motors, 
or the stress-strain diagrams drawn by various types of testing 
machines. The accuracy of the result obtained would in such 
cases be as good as is required, but it is to be remembered that 
graphical methods are not susceptible of great refinement. If, 
therefore, considerable accuracy is required, or indeed in general 
whenever the function is specified by a table of numerical values, the 
method of numerical integration is preferable to the use of the plani- 
raeter. The method of numerical integration must also be resorted to 
when a function is specified by its analytical expression, but cannot 
be integrated in terms of known functions by the methods of 
the Integral Calculus; e.g., the elliptic integrals belong to the class 
of functions which cannot be evaluated by the elementary 
methods of the integral calculus. 

The formulae necessary for numerical integration are derived 
from the formulae already established for interpolation. As in 
the case of interpolation, the order of differences which must be 
taken into account, will depend entirely on the rapidity with 
which the differences decrease as the order increases. It need 
hardly be said that unless the convergence of the series can be 
ensured the process is of no avail. Consider, for example, the 
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determination of the area of a semicircle, the ordinates being 
perpendicular to the diameter. No matter how many ordinates 
are taken (i.e. no matter the order of differences employed), the 
approximation to the area by means of many of the formulae to be 
hereafter established will not be a good one. But if we evaluate 
the area included between an arc of a circle, less than a semicircle, 
a diameter, and the two perpendiculars from the extremities of the 
arc on this diameter, we can, from the same formulae, obtain an 
approximation to the true value with any degree of accuracy we 
please. The reason for the failure in the first instance is that the 
differential coefiicients of the ordinates at the limits of integration 
are infinite, and as a result the series is not convergent. For this 
same reason we must assume that the function represented by the 
observations, and its differential coefficients up to the order of 
differences employed, are continuous ; a condition which is not 
necessarily fulfilled in natural problems, as may be seen by an 
examination of an ordinary barometric curve. 

23. Evaluation of Intep-als by the Integration of 
Infinite Series Term by Term. 

Before applying the methods of finite differences we shall give 
an example showing how in many cases the value of the integral 
of a function which is analytically known may be obtained by 
expanding the function in a uniformly convergent infinite series, 
and then integrating it term by term. 

Example 1. — Consider the complete elliptic integral 



F{k, 



dx 



)o J{\-¥sm'x)' 



f^ 



dx (1 - A'sin^a:)"' 
Jo 

Assuming that | A | < 1 and expanding the term ( 1 - A;^ sin^ x)~^ 

by the binomial theorem, we have 

IT 

Jo *■ 2\2! 

1.3.5,„ . ^ \ 

+ ^r3T^ '''''" + •••/• 
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But rsin-^c^^=l-^-^-(2"-^)'^ 



|"sin=" 
Jo 



2.4.6 ... 2n L'" 
Hence 

^(^, I) = 1(1 + ihr k^ + (^y ^^ + (rlreT^"^ ■••/ ■ 

If A is small this series is very suitable for obtaining the 
value of F ; but it is evident that if k approaches unity a large 
number of terms must be included to get a fair approxima- 
tion to the value of F. 

Let i = 0'L Then to obtain F correct to 4 decimal places it is only neces- 
sary to include two terms of the series. For 

/^Co-1 ') ^i 0-25 0-140625 ) 

'^Voi, 2;= o 1 1 + 100 + 10000 +•■■( 

= ;[{H-00025-H0-0000140625-t-...} 
= :| {1-0025} 
= 1-5747. 

Example ;?.— Find the value of F (k, j) when i = sin 10° 

IT 

Example 3. — Show that I Jl-Jfisin'x is equal to 1-54415 when k has 

the value sin 15°. 

2 fO-6 . 

Example 4.— Show that the value of the integral —pi- I 6-=^ ax is 

V T Jo 

0-60386. 



24. A Formula of Integration Based on Newton's 
Formula. 

We shall now show how integration may be performed by the 
aid of interpolation-formulae. Integrating Ne-wton's interpolation 
formula 
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with respect to n, we have 

I /{a + nw)dn = /„ \dn + 8/. \ndn+ — S^J\ \n{n-l)dn 

+ ^^^n{n-l){n-2)dn + ... 



Hence 



[/{a + nw)dn =/„ + J 8/^ - ^^ S=/, + ^^ S'f^ - ^ S% + ... 
Jo 

^'/{a + nw)dn =/, + 18/,- ^V 8=/, + .V ^Y^ - T^^V 5^. + • • • 
J /(a + nw)dn =f, + i S/j - t\ 8=/ + ^V 8^^ - ^^ ^V. + ■ ■ ■ 

[/(a + nw)dn =/_, + 1 8/_^ - T-V 8^. + ^V 8^,+^ - yVir S^+i + - 

Adding these equations we have 
f{a + nw)dn= /0+/1+ ... +/,_, 

-TV(syi+8y.+ ... +sy,) 
+^(8y„ + 8y^+... +8-y,^,) 

-7VTr(Sy2 + 8y+... +8^+0 

+ 

But 

S/j + 8/3 + . . . + 8/,_, =/ -/„ +/ -/+... +/, -/._, 

=/)■ ~y 

8y: + 8Y,+ ... 8y = 8/^ - 84 + S/. - S/n + ... + 8/^, - S/^.j 
= 8/^+j-6/^ 

sy| + sy^ + . . . + sy,+ J = sy - sy + 8y3 - sy + . . . + sy+, - 8yr 

= sy+.-8y 
8y+sy3+ ... + 8y+, = s»4-8yj+8»yj-8=4+ ... +8%^^-h%^^ 

=6y+^-8yj 



24] NUMERICAL INTEGRATION 67 

Hence 

^J(a + nw)dn = ^f, +/ +/, + . . . +/^_, + y; 

-T¥Tr(sy.+|-sy,)+ 

Putting x = a + nw this gives 

-T¥«(8v;+j-8y3)+ 

This expansion will give the value of the definite integral, pro- 
vided it is possible to obtain values of /(x) outside the limits of 
integration in order to obtain the diiferences. 

When this result is written in the form 



1 /•n+rif 



(x)dx = (/ +/,+ ... +/,_,) + J (./; +/„) - tV (8/;+i - Vl) 
+ #T (8^+1 - 8^1) - r'^ (8^.+^ - sy») + . . . 

it is easily seen that the coefficients are those of x in the expansion 

, ^ : = l+A„x + A^x'' + Ar,!x?+... 

log (1 +05) 

Clausen * has computed the first thirteen of these, viz. : — 

1 , 33953 



= +- 



'~ 12 



il2= + 



24 
^'- 720 



^^"■^160 






3_ 

6C 
863 

60480 
275 

24192 





3628800 


A- 


8183 
■*" 1036800 


A = 


3250433 
479001600 


•"10 = 


4671 
"•" 788480 


^11 = 


13695779093 


2615348736000 


-"12 = 


2224234463 


475517952000 



* Jour, reine angew. Math., 6 (1830), p. 287-9. 
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In the same memoir it is shown that when the limits of 
integration are a- ^w and a + {r + ^)w the coefficients are those 
in the expansion 

J (I +a;)log(l+a:) 

These are also tabulated. 

The formula obtained above must be subjected to another 
transformation when values outside the limits of integration cannot 
be obtained. 

Now 

(i) s/,,+,=8^y;.+8/,_, 
(ii) sy,.+,-2 8Y,+8y,_, = 8^y;. 
i.e. 8=y;.+,-2 8y,= -8y,_, + 8^/; 

(iii) S^fr=S%^. - 3 SVr+i + 3 S'/r-i - Sy^.g 

i.e. 8y,+„ = 8y + 3 8y;^j - 3 8^.., + 8y,_^ 

Substituting these values in the formula we have 

a+rw 

f{x) dx = ^f, +/i +/,+ ... +y,_i + 1/, 

- iV (S/,., - 8/,) - -^ (8y._, + 8y) 

- t¥o (sy.- J - %) - tIo (sy.-. + sy,) 






This form involves only diiferences which can be calculated 
from values of the function which are within the limits of 
integration, and is the form to be used when no values outside 
those limits can be obtained. 

Example 1. — To determine v from the equation 
TT C'^ dx 

using differences of j^. 

1 
Here /^ = 1+^ ! /o=/(0) = l; /,. =/(l) = 0-5 ; 

!0 = 01 ; r = 10. 

Inspection of the formula will show that it is not necessary to form a 
complete table of differences. The highest order of differences retained will 
be the 4th, and as S^f^ and S^'/g are the only two of this order required, the 
central portion of the table may be omitted as in the following scheme ;^- 
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X 


/W 


A 


0-0 


1-00000 


- 990 


0-1 


0-99010 


-2856 


0-2 


0-96154 


-4411 


0-3 


0-91743 


-5536 


0-4 


0-86207 




0-5 


0-80000 




0-6 


0-73529 


-6415 


0-7 


0-67114 


-6138 


0-8 


0-60976 


-5727 


0-9 


0-55249 


-5249 


1-0 


0-50000 




From this table we have 






i/o = 0-50000 






/i = 0-99010 






/a = 0-96154 






/3 = 0-91743 






/, = 0-86207 






/5 = 0-80000 






f, = 0-73529 






/, = 0-67114 






/g = 0-60976 






/g = 0-55249 






4/i„= 25000 






Sum = 7-84982 




lenoe 


1 




• 


19 



A2 

-1866 
-1555 
-1125 



A3 

311 
430 



277 
411 
478 



134 
67 



«4 = 


-0-05249 
-0-00990 


Diff, = 


-0-04259 




0-00478 
-0-01866 


Sum = 


-0-01388 


8Vsi = 


000067 
00311 


Diff. = 


-0 00244 




-0-00067 
0-00119 


Sum = 


0-00U52 



A* 



119 



-67 



-1 (-0-04259) - ^(-0-01388) 



720 



(-0-00244) -j^ (0-00052) 
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= 7-84982 

355 

58 

6 - 1 

= 7-85400 

10 TT 

—^ = 7-85400 
TT = 3-14160. 
Example S. — Determine log. 2 from the equation ^ log, 2 
using differences of ^g. 



13 X dx 
2 1^ 



25. Case when the Upper Limit is not a Tabulated 
Value. 

In the previous section it has been assumed that the quantity r 
in the upper limit was an integer, so that the tabular value of /{x) 
was known when x had the value a + r iv. When the upper limit 
is a + r'w, where r<r' <r+\, and the lower limit coincides with 
one of the values of the argument, so that f{x) is known when 
x = a, a + w, ... a + r IV, a + r+lw, ... the value of the integral 
may be obtained as follows : — 

First determine, as in last section, the value of 

J^ f{x)dx. 

Then, dividing the range (a + ru; a + r'iv) into a convenient 
number of equal intervals, interpolate the values of f{x) for each 
of these sub-intervals. The value of 

Ca+r'w 

f{x) dx 

J a+r 10 

may now be obtained by making use of these new values of /{x) 
and their differences. 

When the lower limit does not coincide with a tabular value of 
the argument, a similar process will enable the value of the integral 
between the given limit and the nearest tabular value of the 
argument to be obtained. 
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Another method which might be employed is to obtain the 
value of I f{x)dx for different ranges of integration, and then 
determine the value for the given range by interpolation. 

(•31° 20' 

Example 1. — Evaluate I ao&xdx, given the values of cosx when 

J 20° 

a: = 20°, 22°, 24°, 26°, 28°, 30°, 32°. 

Forming the table of diflferences we have 



X 

20° 
22° 
24° 
26° 
28° 
30° 
32 



cos X 
0-9396926 

0-9271839 

0-9135455 

0-8987940 

0-8829476 

0-8660254 

0-8480481 



■ 125087 
- 136384 
- 147515 
- 158464 
- 169222 
- 179773 



A2 

-11297 
-11131 
- 10949 
- 10758 
- 10551 

r3o» 



A' 

166 
182 
191 

207 



A« 



16 

9 

16 



First Method.— We shall first of all determine I cos xdx. Substituting 

J 20" 

in the formula 

— {''^"f(x)dx = 4/„ +/,+ ...+/,._, + 4/,. -iii(S/,.-i-5/i)-A(S'-/.-i + SVi) 

'" Ja 



we have 



1 P°° 1 n.^Ros^fiQ a. r -169222-', ,r- 107581 __x!Lr 

___ oosa;da; = 0-4698463-^L^2gog^J ^l_■^_■^a!i^^ ""L- 
0-0349066 j2o» 0-9271839 



1911 
166-1 



P30° 



0-9135455 
0-8987940 
0-8829476 
0-4330127 
= 4-5257896 

cosa:da:=0-1579798. 



120" 



Now by using Newton's formula for backward interpolation 
n( n-\) n(»-l)(»i-2) 

/(a -«i(')=/o -"■*/-!+ 2! *'/-!-' 



3! 



SV-4+- 



determine the values oif{x), i.e. cos x-, when x has the values 30° 20', 30° 40', 
31° 31° 20' 31° 40'. We thus obtain the following table :— 
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X 

30° 0' 
30° 20' 
30° 40' 
31° 0' 
31° 20' 
31" 40' 
32° 0' 



cosa; 
'8660204 

0-8631019 
0-8601491 
0-8571673 
0-8541564 
0-8511167 
0-8480481 



-•29235 

- 295-28 
-29818 
-30109 

- 30397 

- 30686 



A'-i 

-•293 
-290 
-291 

-•288 
-289 



3 
1 
3 

1 



From this table we obtain 



1 f3i°20' f- 30109"! r 

0-^058178j3„ "°"^"" = ;'.8Sim9''^L + 29235j "-L 



27 -A 
r863l019 
0-8601491 
0-8571673 
0-4270782 
= 3-4405189 

fai" 20' 
cos a; cJa; = 0-0200161 
30 



-291' 
-293 



Adding these results together we get 



cos X- cZx = • 1 7 7 99o9 



Second i/ci/iocZ. — Obtain by the general method the values 



"'f 



cos X dx 



when the lower limit is 20° and the upper limits are 22°, •24°, 26°, 28°, 30°, 3^2° 
respectively. We thus obtain the following table : — 



X 


1 cos X dx 


A 


A2 


a? 


20° 


0000000 


325865 






22° 


0-0325865 


321300 


-4565 


-390 


24° 


0-0647165 


316345 


-4955 


-385 


20° 


0-0963510 


311005 


-5340 


-381 


28° 


0-1274515 


305284 


-5721 


-370 


30° 


0-1579799 


299193 


-6091 




32° 


0-1878992 
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Using Newton's formula of backward interpolation we obtain the value of 
the required integral, viz. 

f3i<>20' i(j-n ■i,l),-^\(^-'2\ 

oosa;rfx- = 0-1878992-J(299I93) + ^^4rr-(-6091)-^^^-4T^ ^(-370) 

J 20" ■" ■ *5 ! 

= 01779961. 

The values obtained by the two melhods agree to 6 places of decimals. 
The result correct to 8 places is 017799598, but as the last figure is always 
forced, the one method can hardly be assumed better than the other. Probably 
the latter method is to be preferred, as in the former the division of the range 
between the given upper limit and the new upper limit into suitable intervals 
may give rise to values of the argument consisting of many decimal figures 
which would increase the amount of labour involved in the interpolation. 

C 30" 20' 

Examples. — Evaluate I sin.rd.t, given the values of sina; when 

J 20» 

a; = 20°, 22°, 24°, 26°, 28°, 30°, 32°. 

Ci-o dx 1 

Examples. — Evaluate I . , if the values of ,q are given for 

integral values of x. 

26. A Formula Based on Bessel's Expansion. 

Another series for the definite integral is obtained by integrating 
Bessel's formula of interpolation, 

, , ,< , 11 in- ]) -„ , n (n - I) (n - i) 
/(a + n lu) = f.f, + {n - i) Sf, + ^^, M S% + -^— 3T ^ ^ A 

n { n- - I) {71 - 2} 
+ j^ ■/xSy,+ ... 

Putting x = a + nw we have by integration 

1 f"+"' f' 

— I f\x)dx= I f{a + nw)dn 

10 Ja ' Jo 

fi n o, .Pn(ti- 1) , 

= fx/A dn + 8/ {n-)^dn + ,xS-jA ^, dn 
-in •'0 Jo ^ • 



■"O 



■"'•^j. 



■1 n {n - 1)(m-A) 
3! 



4! 
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Similarly, 

— f{x)dx ^ /x/^ - ^V M ^y* + -Nts /* Sy^ + . . . 

- f{x)dx = i..y;._j-3-V/.sy,,..+^Vu/^5yr-^+... 

*" J «+()■— 1)10 

Hence 

-J ^^ /(x) c/ 33 = /.y^ + ^4 + . . . + /./^_^ 

+ tVV { i^-^V^ + A. sy^ + . . . + /. sy _ J 
+ 

- t\ { i 8y„ + sy + 8y + . . . + sy_, + ^ sy } 
+ tVV {i sy„ + sy + sv; + . . . + sy,_, + 1 sy,} 
+ 

or, since 

sy + sy+.-. + sy-: = s/„-5y+84-sy„ + ... + 8/,,_^-S/^_„ 

and therefore 

sy + sy+... + sy_ = 8'/,..,-% 

we have 

■J f{x)dx = i/o +/i +/2 + . . . +/._! + \fr 

-TV{|syo-sy + s/,., + i8y,} 
+ Tytr{i8Vo-8y,+8=y_^+|8y} 
+ 

But 

iSyo-8/^=|8/^-J8/_j-6/^ 



w 
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and, similarly, 

Hence, finally, we obtain the formula 

— 1^ fix) dx = \f, +/, +/,+ ... +/._, + y; 

- t\ {f^ s/, - /. &/,) + J^ {,, sy, - /. sv„) 
+ 

Example i.— Calculate I — to 8 places of decimals. 
J 100 ^ 

Here f{x) = — , /„=100, r = 5, w = \. 

In order to obtain the differences required, it will be necessary to form 
f(jo) for all integral values of x from a; = 98 to x = 107. These may be 
obtained from Barlow's Tables. We thus obtain the following :- 



X 

98 
99 
100 
101 
102 
103 
104 
105 
106 



0010204082 
0010101010 
0-010000000 
009900990 
0-009803922 
0-009708738 
0-009615385 
0-009523810 
009433962 



A 

- 103072 

-101010 
( - 100010) 

- 99010 

- 97068 

- 95184 

- 93353 

- 91575 
(- 90712) 

89848 

- 88168 



A2 

206-2 
2000 
1942 
1884 
1831 
1778 
1727 
1680 



A' 



-62 

(-60) 
-58 

-58 

-53 

-53 

-51 

(-49) 

-47 



107 0-009345794 

Substituting in the formula 

105 dx 

X 



J 100 
we get 



i /« + .A +A + L+J\ + U,-h{y. SA - M «/o) + -Vo (/x sVs - M «Vo) 



(•105 dx 

1 — = 048790940 - -000000775 + -000000000 

J 100 ^" 



= 0-018790165. 
The value correct to 10 places is 0-0487901642. 
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Note. — It is to be remarked that the expression 

i/o+/i +/2+/3+/4 + U=0-048790940 

gives the answer correct to six decimal places ; this approximation is quite 
sufficient for most practical purposes. 

f2 dx 

Example 2. — Approximate to the value of the integral I 113? -'23?+ 2 ) 

27- The Euler-Maclaurin Expansion. 

In the foregoing the value of the integral is expressed as a 
series involving differences. In many cases in which the function 
is known analytically, it will be found more convenient to use what 
is known as the Euler-Maclaurin expansion which employs 
differential coefficients instead of differences. The rapidity of 
convergence is of the same order as that given by Bessel's formula. 
It may be obtained as follows ; — 

Since 

= ^y{a + w)-f[a-w)} 

= n/W + «'/'(«)+^/" («)+■■■ 

W~ 

-f{a) + wf{a)- ;^/" («) + ...} 



--wj («'^+ 3I/ (a)+... 



and similarly 






IX 87, = w^f" {a + rw)+ — /"' (a + riu) + 



we have by substituting in the formula 

1 fl a+r if 






fix) dx=y, +y, +/„+... +/._! + y. 
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the result 

^\7' " -^ ^"^ «^* = i /o +/i + • • ■ +/.-: + J/. 

- tV ^ ifr' -/o') + rk ^' (/■'" -/o'") 

This is the EulerMaclaurin* expansion. The coefficients of the 
various items in it are proportional to the well-known Bernoullian 
numbers. For iff we put /(x) = e"^-" in the formula and at the 

same time change the limits to a - w and a + w we get 

\ ro+io 

— e'- dx = i/„ +/, + hA + A w ifJ -//) + Bir (//" -/„'") + . . . 

"-' J a— to 

where A, B, ... are constants. 
But 

— I e*~" dx '=^ — I e' dx 

W J<i_,i, W J -10 

= — le" - e-'°). 
w 

Hence 

— (e" -«-") = i/„ +/, + iA + A w {/; -/„') + B v^ (//" -/o'") + ■ • • 

= I e-" +\+\e^ + Aw (e" - «-") + B ii? (e" - «-'") + . . . 
from which it is easy to deduce that 



•2 «""■' - e-'""- 



1 -Aw'-B^d" 



— cot — . = I - Aw^- Bw* - ... 
2i 2t 

But the Bernoullian numbers are the quantities 5„ B.,, ... in the 
expansion 

-r^ 2^ ?j6 

and therefore i, B, ... are proportional to B„ B^, ... 



* For a rigorous discussion, c/. Whittaker and Watson's Modern 
Analysis, p. 128. 

t Poisson, M^m. de I'Acad. des Sc, li523. 
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Example i.— Calculate I "^ t° § places of decimals, and check the 

result by computing it as -\oge{.\-^^) by the logarithmic series. 

Here /o=^, «"=!, '" = 5- 
Then 

Jj5 "^"T'95 96 97 98 99 2 100 



T'^V 1002"^95V"*"i-" V 100''"^95V 



1 / - 0-O00100000\ , 1 / - 0000000010\ 
0-005263158 - t -1(^ + -0001 10803 i "'^V +0-000000012^ 



010416667 
0010309278 
0-010204082 
0-010101010 
0-005000000 

: 0-051294195 
-0-000000900 

: 0-051293295 



Again 



(0-05)^ (0-05)^ (0 05)^ (0-05)° (0-05)° 



= 0-050000000 

0-001250000 

0-000041667 

0-000001563 

0-000000063 

0-000000003 

= 0-051293296 

The two results will be seen to agree to 8 places of decimals. 

f 105 dx 
Example 2. — Evaluate I — to 8 places of decimals. 

Iioo ^ 



Example 3. — Approximate to the value 



28. An Exceptional Case. 



of 1 cos- X dx. 
J 



Legendre* remarked that if the odd differential ooeiBcients above b. 
certain order take the same value at both limits, the formula fails to give an 
accurate value of the definite integral. For instance, if the integral under 

* Fonctions Elliptiques, Vol. 11., p. 57. 
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discussion were the elliptic integral of the secood kind taken between limits 
and ir/2, viz. 

I ij 1— A-^sin^ xdx, 
Jo 

all the odd differential coefficients become zero at both limits, and thus the 
value of the integral would appear to be given by 

But the value given by this expression is not in close agreement with the 
true value of the integral. The reason of this is that the numerical 
coefficients introduced in the successive differentiations increase without 
limit, 80 that each term takes the form oo x 0, which is, of course, inde- 
terminate. 

Better methods of evaluating the three kinds of elliptic integrals will be 
given in another tract in the present series. 

29. The Formulae of Woolhouse and Lubbock. 

We shall next derive from the Euler-Maclaurin expansion a 
formula which is of use when it is required to evaluate the sum of 
a large number of terms. 

If the interval between consecutive values of the argument in 
the Euler-Maclaurin formula is w, we have 

- f "'"' /H dx = /„+/ +/,+ ...+/,- i (/„+/) - i^w (/; -/»') 

to Ja 

+ A-T!^(fr"'-/o")+- 

and if this interval is, w f m the formula becomes 

- f "^"'/(«) rf a; =/„ +/i +/ ^ +...+/.- i (/o +/.) 
rv J a '" '» 

- A ~ ifr -fo') + TW ~3 if'-'" ~fo"') + ■■■ 

Subtracting the latter from m times the former we deduce 

_ 1 

/o +fl +/l +■••+/■ = »« (/o +/l + • • • +fr) - —2- (^0 +/^) 

~ 12 m^-^' ^'^ 
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This is Woolhouse's Formula. Obviously by using the right- 
hand expression instead of the left-hand one, the labour involved in 
determining the sum/o-|-/_i^ +/ 2 + ... +fr is enormously reduced. 

A somewhat similar formula, in which differences are employed 
instead of diiferential coefficients, is due to Lubbock, Expressing 
the derived functions in Woolhouse's formula in terms of the 
successive differences of /j, J\, ... fr we obtain the formula in 
question, 

/o +/ 1 +/i + • • • +/. = m (/o +/i + • ■ • +./;.) - ^^ (/o +/.) 

(m°-l)(19>."--l) 

720 m^ ^ -^'-t ^JV 



(m°--l)(9TO^-l) 
480 w' 



(s^z-.+sy^) 



350 1 

Example 1. — To determine 2 — • 

»l = 30O " 

Let m=10, w) = 10. 

Then from Woolhouse's formula we have 

^„„ n =^0\300 310 320 330 340 350/ 

n = 300 ' 

2 \300 350/ 12 V 350'^ 300V 

= 0-18512761 -0-02785714 
2432 
= 0-15724615 

Using Lubbock's formula we get 

360 1 
)t = 300 



_ 99 . 1899 fS3/-.. 
720000 ^ •'k 


-«»/») 


= 0-18512761- 


002785714 






1938 






487 






3 






3 


= 015724616 







. 9 » ■ 3 99 
480000 



W,+^%) 
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To test the accuracy of these answers, the reciprocals of the numbers 

300, 301, 302 350 were taken from Barlow's tables, and summed by 

means of a comptometer. The result was found to be 0-15724616. The 
closeness of the approximation is remarkable. 

175 1 
jErampie g. — Obtain the value of V ^• 

n=100 " 

Example ?,. — By giving n increasing positive values, show that 



L 



1 

(l + i + 4+ ••• + — -log„?i) = 0-57721566.. 



30. Formulae of Approximate Quadrature. 

The expansions which have been obtained in the preceding 
articles furnish the value of an integral to any required degree of 
accuracy, provided a sufficient number of terras is taken. We shall 
now consider certain formulae which are in their nature only 
approximate, and which, though frequently useful, cannot be used 
where great accuracy is required. 

Suppose that it is required to integrate a function f(x) between 
limits, which may, without loss of generality, be taken to be - 1 
and + 1. We shall endeavour to obtain expressions of the type 

H,f{K) + HJ(S) + . . . + //„/(/i„) 

which represent the integral as closely as possible, where \, A„ ... A„ 
are (n + 1) values of x within the range of integration and Ao, A,, ... A„, 
Z?oi ii\i ••■ Hfi do not depend on the function/(a;). 

Let us first choose //„, H-^, ... H,„ so as to make the formula 
strictly accurate, so long as f{x) is a polynomial of degree less 
than n + l. 

Let F(x) = {x- ho) {x - A,) ... (x- h„). 

Fix) 
Then — —^ is a polynomial of degree less than n+l, so that the 

x-h,. 

Fix) 
formula must be strictly accurate when for f{x) we put — . 

This gives 



-1 X- 



= H,F'{K). 
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Hence, whatever h„, h^, . . . h„ may be, in order that the above con- 
dition may be satisfied, the values of H^, IIj, ... ff„ must be 
given by 



F'{h,.) j.ix-/,. 



In practice it is convenient to make the intervals between 
successive values of the argument equal to each other. Suppose, 
then, that /t„, h^, ... h„ are chosen so as to divide the whole range 
of integration into n equal parts. 

Then A„ = - 1 

h, = -1+ — 
n 

4 

/io = - 1 + — 

n 



H = 



[h,. - Ao) . . . {K - A,._,) {h,. - A,+, )...{K- K) 

C+1 

X I (a; - Ao) ... (a; - A,._,) {x - h,.^^) ... (x- h„)dx 

1 

^2r 2 (r - 1) _2_ -2^ ^1 -2(n-r) 
n n n n 71 n 

( - 1)"-' 



Now, let t=-~ {x-\- 1). 
Then 

H,= -^-rj r-. t{t~\)...{t-r+\){t-r-\) ...{t-n)dt. 

Hence, finally, if we put 

hQ = a, ... h^ = a + rw 
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we obtain 



1" 



fo+nw H 

f{x)dx = ^ H,f{a + rw) 

1=0 

where 

"'^ :;Ti~;:^,XHt-\)...(t-r+\){t-r-\)...{t-n)dt. 

This is known as Cotes' Formula of Quadrature. By means of 
it we maj' calculate the value of the integral without first of all 
forming a table of differences. 

Special Cases. 

(i.) Letn=l. Then 

\.f(x) dx = HJ(a) + H,/(a + w) 

J a 

where ir„ = ^-j-yl (t-l)dt = lw 

fo+w -J, 

f{x)dx = -{f{a)+f{a + w)}. 



Hence 



This is known as the Trapezoidal Rule. 
(ii.) Let»i = 2. Then 

[^' "fix) dx = If of {a) + HJ{a + w) + H,f{a + 2w) 



where 

0! 2 

n 

t(t--2)dt =|io 



Ho=K^,\\i--^){t-2)dt = ^w 



^'^iTT-ilo 



^==2-!^!'" '^'-'^'^' 



So '^' 



Hence 

/(a;) (fee = -/(a) + —f{a + w) + -f(a + 2w). 

This is called Simpson's * Bule or the Parabolic Bute. 



It was, however, first given by James Gregory. 
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Corollary. — If we divide the whole range of integration, say 
(a, 6), into 2 n equal parts, and apply Simpsons Rule to each of 
them, the value of the integral may be expressed as 

I = -g^ {2/0 + y^,, + 2 (2/0 + yj + . . . + 2/2,._;) + 4(2/1 + 2/3+...+ 2/.2„_i) } 
where y„, 2/1, ... 2/2,1 are the values of /(as) when x has the values 
a, a + — — , ... 6. This is Simpson's Rule in its generalised 

form. The rule may be stated thus ; — 

Divide the area into an even number of strips by equidistant 
ordinates 2/01 J/u ••• 2/2" • ^^ '^* *'*''* ^f '^^ extreme ordinates add 
twice the sum of the other ordinates with even suffixes, and four 
times the sum of the ordinates with odd suffixes, and multiply this 
total by one-third of the commoti distance between the ordinates* 

(iii) Let n = 3. Then 

fa+Sw 

f{x) dx = HJ{a) + HJ{a + w) + HJ{a + 2w) + HJ{a + 3 w) 

J a 

where ^0 = oyT! Jo ^^ ~ ^^ ^* ~ ^^ ^*' ^^ '^'^^"' 

= 271111 '(«-!) (''-3) <^« -l^ 



H, 



H,= ^ '" 



Hence 

V^'°f(x)dx = |t«{/(o) + 3/(ci + «)) + 3/(a + 2w) + /(a + 3M))}, 

which is often termed the Three-Eighths Rule. 

(iv) Let « = 6. Then it may be shown, as in the above 
special cases, that 

(•0+6 w 

/(ic) c^!B = T-i^ {41/(a) + 216/(a + M>) + 27/(a + 2 w) 
•'" + 272/(a + 3 w) + 27/(a + 4 m;) + 216/(a + 5 w) 

+ 41/(a + 6w)}. 

* The reader is warned that the notation varies in the different text- 
books, so that in all cases an examination of the notation used should first 
of all be made, 
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Adding 

ths ^/i<^ + 3 w) = xk {/(») - 6/(a + «>) + 1 5/(a + 2w) 

-■20/(a+3w) + l5f{a + iw)-6/{a + 5w) 
+f{a + Giu)} 

we obtain 

r "/(^)'^« + TTirSV(« + 3w) = xTTr{W(«) + 210/(a + M)) 

+ 42/(a + 2 to) + 2o2/(a + 3 w) 
+ 42/(a + iw) + 210/(a + 5 w) 
+ i2/{a + Qw)} 

= A {/(«) + 5/(a + w) +/(« + 2 w) 
+ 6/(a + 3 w) +/(o 4- 4 M)) 
+ 5/((i + 5w)4-/(o + 6 w)}. 

Neglecting the term -[-|^ 8°/(a + 3 w), which will be fairly 
small, we get Weddle's Rule of Quadrature for 6 intervals. 

y (a;) rf a; = -i\ { /(a) + 5/(a + 1«) +/(a + 2 to) + 6/(a + 3 to) 
+/(a + 4 zo) + 5/(a + 5 w)+/(a + 6 w) }. 

31. Comparison of the Accuracy of these Special 
Formulae. 

We may compare the accuracy of these formulae by evaluating, by each 

f7 dx 

of these methods, the definite integral I — , the correct value of which is 

Jl ^ 
logj = 1-94591. 

(i) Using the Trapezoidal Rule six times, i.e., dividing the region into 
six equal portions and applying the rule to each of them, we get 

[ -^ = 4/((i)+/(a + i(;)+/(o + 2«))+ ... +/(a + 5«)) + i/(a + 6i«) 

= i-l+Hi + i + ^ + i-i 
= 2-02143. 



(ii) Using the Parabolic Rule repeated three times we get 

f - 



'' — = i/(o) + |/(a + iu)+§/(o + 2w)) + */(a + 3M>) + »/Ca + 4t«) 



= ^(l+l + l + T + f + F + T) 

= 1-95873. 



86 INTERPOLATION AND NUMERICAL INTEGRATION [ch. iv 

(iii) Applying the Three-Eighths Rule twice we have 

[' -^ =a{/(a) + 3/(a + .o) + 3/(a-H2«,) + 2/(a-l-3u>) + 3/(a4-4«7) 
J^ " +Zf[a + 5w)+f(a + Qw)) 

= |{l+t + l + f + l + T + T} 
= 1-96607. 

(iv) Weddle's Rule gives 

= 1-95286. 

None of these gives the result with any high degree of accuracy, the 
errors being 

(i.) 0-07552 
(ii.) 0-0128-2 
(iii.) 0-02016 
(iv.) 0-00695 

n dx 

Example 1. — Evaluate I ~~rk hy each of the methods of § 30. 

Example 2. — Plot the curve j'^ = 2a: + 25 and find the area in square inches 
between the curve, the j/-axis, and the double ordinate at a;=12, the unit 
being one inch. 

Example 3. — The ordinates of the boundary of the deck of a, ship are 
8, 30, 36, 40, 42, 42, 40, 38, 34, 28, and 8 feet respectively, and the common 
distance between them is 21 feet. Find the area of the deck. 

Example ^. — The velocity of a train which starts from rest is given by the 
following table, the time being reckoned in minutes from the start and the 
speed in miles per hour. 

t 2 4 6 8 10 12 14 16 18 20 

mlh 10 18 25 29 32 20 11 5 2 at rest. 
Estimate the total distance run in 20 minutes. 



32. Gauss' Method of Quadrature. 

The accuracy of the method given in § 30 may be inci-eased if 
we no longer suppose the intervals between the successive values 
/t„, /ii, ... A„ of the argument to be equal. 

Let us now endeavour to choose H^, II^, ... II, ^, //„, A,, ... h„ in 
such a way that the formula 

V^ f{x)dx = HJ[K) + HJ{K)+ ... +HJ\K) 
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shall be accurate for all functions f{x) of degree equal to or less 
than (2»i+l). 
Taking 

f{x) = F{x)4>{x) 

where <^ (a;) denotes any polynomial of degree <(«.+ 1), and where 

F (x) = {x- hfi) {x - hy) (x-L) ... {x - A„) 
we see that the condition to be satisfied may be written 

f F{x)^{x)dx = 0. 

Now it is well known that if P,,^., {x) denotes the Legendre 
polynomial* of degree (n+ 1), we have 

r Pu+^{x)ci^{x)dx= 

and, conversely, this condition is sufficient to determine the 
Legendre polynomials. Hence we have (save for a multiplicative 
constant) 

F{x) = P„+,{x) 

and therefore the quantities Jiq, /tj, ... A„ must be the roots of the 
equation 

F„+Ax)^0. 

The following are the values of the constants h^, h^, ... h^, H^, Hi, ... H ^^ 
in the simplest cases. 



^,= 1 



(i) When 71 = 1. 




\ = 


1 

-V3 


K = 


1 

v'3 


(ii) When Ji = 2. 




;,„ = 


- -Jl 


'-1 = 





K = 


VI 






* Of. Whittaker and Watson's Modem Analysis or Byerly's Fourier's 
Series and Spherical Harmonics. 
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(iii) Whenjj = 3. 

''1 = - V(? - -i^ n/30) H, = ^ + ^ V30 

h = >/(?-- A V30) ^2 = H sV x'30 

h-i = V(f + 3F V30) H-, = h-i^ V30 

(iv) When 7i=4. 

322- 13 ^'70 

K = v'(-5 + 1% V70) ^„ = -^^"--- 

128 
^2 = H^ = 225 

/- " Mn> rr 322+13s/70 

322-13\/70 



Example 1. — Determine the value of I qT"*- 
(i) Using 3 ordinates 



900 



1 



'+^ _d«_ ^ _5_ 1 ^ J_ + A 1 

l3+< "9""3- s/3 9 • T 9 '3+^/5 



1 



= 0-69312165. 
(ii) Using 5 ordinates 

'+^ dt _ 322-13s/70 1 322+13^70 1 

3+<~ 900 ■ 3 - ^/(| + e^a v'TO) 900 's-^d-^^yo) 

128 J_ 322+13V70 1 

+ 225' 3 + 900 3+ V(«-AV70) 

322-13V70 . I 

+ 900 3+^(1 + ^^70) 

= 0-693147157. 

The correct value ia 0-69314718, so that using 3 ordinates the error is in the 
5th decimal place : using 5 ordinates, in the 8th decimal place. 



f+l dx 
Example 2. — Determine the values of (i) I TTo ' 

f+l dx 
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MISCELLANEOUS EXAMPLES. 

1. The values of a function /(a;) for the values 1050, 1-060, 1-070, 1080, 
1-090, 1-100 of the argument are 1-25386, 1-26996, 1-28619, 1-30254, 1-31903, 
1*33565. Use Newton's formula of quadrature to shovr that the value of 

f 1-100 

/(a;)da- 
Jl-OSO 

is 0-06472. 

2. The ordinates of a plane curve are of the following lengths, 3-5, 4-7, 
5-8, 6-8, 7-6, 8-1, 8-0, 7-2 feet, and they are 4 feet apart. What is the area 
included between the two end ordinates? and what is the distance of the 
centre of gravity of that area (i) from the extreme left-hand ordinate, and 
(ii) from the base line ? 

3. Approximate to the values of the following definite integrals ; — 

(i) I dxsJ{Z-x + x^) 

cos X 
(ii) I dx 



i 



dx 



Va-^sin^a;)" 

4. Use (i) Lubbock's formula, (ii) Woolhonse's formula, to determine the 
value of an annuity-certain for 36 years, interest being at the rate of 3%. 

5. Verify that the area of the curve y = A + Bx+C3fi + I>3^ between the 
limits x=A and x= -h is equal to the product of h into the sum of the 
ordinates at a; = A/v'3 and x= -ft/^/3. 

In the case of the curve y=f{x) = A + Bx + Cx' + Dx? + Ex' + Fx^ verify in 
like manner that the area between a: = A and a:= - A is equal to 

{5f{hJl)+S/(0) + 5f{-hJj)}li/9. 

6. Corresponding values of x and y are given iu the following table, the 
unit being the inch. Find the volume of the solid of revolution formed when 
this curve revolves about the a;-axis. Find also the centre of gravity of this 
uniform solid. 

^ 4 8 12 16 20 

y 0-50 2-58 4-41 5-40 5-10 050 
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7. A body of 200 lbs. is acted on by a variable force F. No other force 
acts on the body. When the body has passed through the distance x feet 
the force in pounds is as follows ; — 



X 


0-1 


0-2 


0-3 


0-4 


0-5 


0-6 


7 


0-8 


0-9 


1 


F 20 


21 


21 


20 


19 


18 5 


180 


13-5 


9 


4-5 






Determine the work done upon the body from x = Q to a; = 0'4. If the 
velocity is when a; = 0, what is its value when a;=0"4 ? 

8. The semi-ordinates in feet of the deck plan of a ship are respectively 
3, 16-(i, 25-5, 28-G, 298, 30, 29-8, 29-5, 285, 242, and 6-8, the common 
interval between them being 28 feet. Find the area of the deck. 

9. The "half" ordinates in feet of the midship section of a vessel are 
12-5, 12-8, 12-9, 13, 13, 12-8, 12-4, 11'8, 104, 6-8, 0-5 respectively, and the 
common interval is 2 feet. Find the area of the whole section and the 
position of its centre of gravity. 

10. The area of a ship's load-water plane is 8000 sq. ft. ; the body 
below the load-water plane is divided into six portions by equidistant 
horizontal sections, 3 feet apart, whose areas in sq. ft. are respectively 
7600, 7000, 6000, 4500, 2800, and 250. Find (a) the displacement in tons, 
and (6) the number of tons which must be taken out of the ship to lighten 
her 4 inches. 
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